summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
Diffstat
-rw-r--r--CMakeLists.txt2
-rw-r--r--src/CMakeLists.txt12
-rw-r--r--src/kernels/CMakeLists.txt45
-rw-r--r--src/kernels/templates/center_template.in14
-rw-r--r--src/kernels/templates/common.in38
-rw-r--r--src/kernels/templates/definitions.in6
-rw-r--r--src/kernels/templates/lamino_template.in14
-rw-r--r--src/kernels/templates/roll_template.in14
-rw-r--r--src/kernels/templates/z_template.in19
-rw-r--r--src/kernels/tools/make_burst_kernels.py63
-rw-r--r--src/lamino-roi.c153
-rw-r--r--src/lamino-roi.h23
-rw-r--r--src/ufo-lamino-backproject-task.c828
-rw-r--r--src/ufo-lamino-backproject-task.h53
14 files changed, 1280 insertions, 4 deletions
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 54e17ec..19602f2 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -21,6 +21,8 @@ include(PkgConfigVars)
configure_paths(UFO_FILTERS)
set(PKG_UFO_CORE_MIN_REQUIRED "0.6")
+# Backprojection burst mode, must be one of 1, 2, 4, 8, 16
+set(BP_BURST "16" CACHE STRING "Number of projections processed in one pass")
option(WITH_PROFILING "Enable profiling" OFF)
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index bfb370f..c4cd2d5 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -34,6 +34,7 @@ set(ufofilter_SRCS
ufo-get-dup-circ-task.c
ufo-ifft-task.c
ufo-interpolate-task.c
+ ufo-lamino-backproject-task.c
ufo-loop-task.c
ufo-map-slice-task.c
ufo-measure-sharpness-task.c
@@ -94,6 +95,9 @@ set(ifft_misc_SRCS
set(retrieve_phase_misc_SRCS
common/ufo-fft.c)
+set(laminoaux_SRCS
+ lamino-roi.c)
+
file(GLOB ufofilter_KERNELS "kernels/*.cl")
#}}}
#{{{ Variables
@@ -107,6 +111,7 @@ if (CMAKE_COMPILER_IS_GNUCC OR ("${CMAKE_C_COMPILER_ID}" STREQUAL "Clang"))
"-Wno-missing-field-initializers -Wpointer-arith "
"-Wredundant-decls -Wshadow -Wstrict-prototypes -Wwrite-strings")
endif()
+add_definitions ("-DBURST=${BP_BURST}")
#}}}
#{{{ Dependency checks
find_package(TIFF)
@@ -222,9 +227,8 @@ foreach(_src ${ufofilter_SRCS})
ARCHIVE DESTINATION ${UFO_PLUGINDIR}
LIBRARY DESTINATION ${UFO_PLUGINDIR})
endforeach()
+#}}}
-# copy kernels
-foreach(_kernel ${ufofilter_KERNELS})
- install(FILES ${_kernel} DESTINATION ${UFO_KERNELDIR})
-endforeach()
+#{{{ Subdirectories
+add_subdirectory(kernels)
#}}}
diff --git a/src/kernels/CMakeLists.txt b/src/kernels/CMakeLists.txt
new file mode 100644
index 0000000..0802482
--- /dev/null
+++ b/src/kernels/CMakeLists.txt
@@ -0,0 +1,45 @@
+cmake_minimum_required(VERSION 2.6)
+
+# make burst backprojection kernels
+set(COMMON_TEMPLATE templates/common.in)
+set(DEFINITIONS templates/definitions.in)
+set(Z_TEMPLATE templates/z_template.in)
+set(LAMINO_TEMPLATE templates/lamino_template.in)
+set(CENTER_TEMPLATE templates/center_template.in)
+set(ROLL_TEMPLATE templates/roll_template.in)
+set(GENERATOR tools/make_burst_kernels.py)
+set(Z_KERNEL z_kernel.cl)
+set(LAMINO_KERNEL lamino_kernel.cl)
+set(CENTER_KERNEL center_kernel.cl)
+set(ROLL_KERNEL roll_kernel.cl)
+
+add_custom_command(
+ OUTPUT ${Z_KERNEL} ${CENTER_KERNEL} ${LAMINO_KERNEL} ${ROLL_KERNEL}
+ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/${GENERATOR} ${GENERATOR}
+ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/${COMMON_TEMPLATE} ${COMMON_TEMPLATE}
+ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/${DEFINITIONS} ${DEFINITIONS}
+ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/${Z_TEMPLATE} ${Z_TEMPLATE}
+ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/${CENTER_TEMPLATE} ${CENTER_TEMPLATE}
+ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/${LAMINO_TEMPLATE} ${LAMINO_TEMPLATE}
+ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/${ROLL_TEMPLATE} ${ROLL_TEMPLATE}
+ COMMAND python ${GENERATOR} ${Z_TEMPLATE} 1 2 4 8 16 > ${Z_KERNEL}
+ COMMAND python ${GENERATOR} ${CENTER_TEMPLATE} 1 2 4 8 16 > ${CENTER_KERNEL}
+ COMMAND python ${GENERATOR} ${LAMINO_TEMPLATE} 1 2 4 8 16 > ${LAMINO_KERNEL}
+ COMMAND python ${GENERATOR} ${ROLL_TEMPLATE} 1 2 4 8 16 > ${ROLL_KERNEL}
+ DEPENDS ${GENERATOR} ${COMMON_TEMPLATE} ${DEFINITIONS} ${Z_TEMPLATE} ${CENTER_TEMPLATE} ${LAMINO_TEMPLATE} ${ROLL_TEMPLATE}
+ COMMENT "Generating burst backprojection kernels")
+
+add_custom_target(burst ALL
+ DEPENDS ${Z_KERNEL} ${CENTER_KERNEL} ${LAMINO_KERNEL} ${ROLL_KERNEL})
+
+install(FILES ${CMAKE_CURRENT_BINARY_DIR}/${Z_KERNEL} DESTINATION ${UFO_KERNELDIR})
+install(FILES ${CMAKE_CURRENT_BINARY_DIR}/${CENTER_KERNEL} DESTINATION ${UFO_KERNELDIR})
+install(FILES ${CMAKE_CURRENT_BINARY_DIR}/${LAMINO_KERNEL} DESTINATION ${UFO_KERNELDIR})
+install(FILES ${CMAKE_CURRENT_BINARY_DIR}/${ROLL_KERNEL} DESTINATION ${UFO_KERNELDIR})
+
+# copy kernels
+file(GLOB ufofilter_KERNELS "*.cl")
+
+foreach(_kernel ${ufofilter_KERNELS})
+ install(FILES ${_kernel} DESTINATION ${UFO_KERNELDIR})
+endforeach()
diff --git a/src/kernels/templates/center_template.in b/src/kernels/templates/center_template.in
new file mode 100644
index 0000000..101b869
--- /dev/null
+++ b/src/kernels/templates/center_template.in
@@ -0,0 +1,14 @@
+ pixel.x = mad(voxel.x, cosines{1}, mad(voxel.y, sines{1}, x_center_current));
+ pixel.y = mad(tmp_x, sines{1}, mad(tmp_y, cosines{1}, tmp));
+ rotate ();
+ result {2}= read_imagef (projection_{0}, sampler, pixel).x;
+%nl
+ float x_center_current = mad((float) idz, x_center.y, x_center.x);
+ tmp = 0.5f + ((int) x_center_current) - x_center_current;
+ voxel.x = mad((float) idx, x_region.y, x_region.x + tmp);
+ voxel.y = mad((float) idy, y_region.y, y_region.x + tmp);
+ tmp = mad(z_region.x + 0.5f + ((int) y_center) - y_center, sin_lamino, y_center);
+ tmp_x = voxel.x * cos_lamino;
+ tmp_y = -voxel.y * cos_lamino;
+
+{}
diff --git a/src/kernels/templates/common.in b/src/kernels/templates/common.in
new file mode 100644
index 0000000..f31ab64
--- /dev/null
+++ b/src/kernels/templates/common.in
@@ -0,0 +1,38 @@
+kernel void backproject_burst_{0} (
+{1}
+ global float *volume,
+ const sampler_t sampler,
+ const int3 real_size,
+ const float2 x_center,
+ const float y_center,
+ const float2 x_region,
+ const float2 y_region,
+ const float2 z_region,
+ const float2 lamino_region,
+ const float2 roll_region,
+ float sin_lamino,
+ float cos_lamino,
+ const float{2} sines,
+ const float{2} cosines,
+ const float norm_factor,
+ float sin_roll,
+ float cos_roll,
+ const int cumulate)
+{{
+ int idx = get_global_id (0);
+ int idy = get_global_id (1);
+ int idz = get_global_id (2);
+ float result, tmp, tmp_x, tmp_y;
+ float2 pixel;
+ float3 voxel;
+
+ if (idx < real_size.x && idy < real_size.y && idz < real_size.z) {{
+{3}
+
+ if (cumulate) {{
+ volume[idz * real_size.x * real_size.y + idy * real_size.x + idx] += result * norm_factor;
+ }} else {{
+ volume[idz * real_size.x * real_size.y + idy * real_size.x + idx] = result * norm_factor;
+ }}
+ }}
+}}
diff --git a/src/kernels/templates/definitions.in b/src/kernels/templates/definitions.in
new file mode 100644
index 0000000..4d71450
--- /dev/null
+++ b/src/kernels/templates/definitions.in
@@ -0,0 +1,6 @@
+#define rotate() pixel.x -= x_center.x; \
+ pixel.y -= y_center; \
+ pixel.x = pixel.x * cos_roll + pixel.y * sin_roll; \
+ pixel.y = -pixel.x * sin_roll + pixel.y * cos_roll; \
+ pixel.x += x_center.x; \
+ pixel.y += y_center;
diff --git a/src/kernels/templates/lamino_template.in b/src/kernels/templates/lamino_template.in
new file mode 100644
index 0000000..e3453b3
--- /dev/null
+++ b/src/kernels/templates/lamino_template.in
@@ -0,0 +1,14 @@
+ pixel.x = mad(voxel.x, cosines{1}, mad(voxel.y, sines{1}, x_center.x));
+ pixel.y = mad(tmp_x, sines{1}, mad(tmp_y, cosines{1}, tmp));
+ rotate ();
+ result {2}= read_imagef (projection_{0}, sampler, pixel).x;
+%nl
+ tmp = 0.5f + ((int) x_center.x) - x_center.x;
+ voxel.x = mad((float) idx, x_region.y, x_region.x + tmp);
+ voxel.y = mad((float) idy, y_region.y, y_region.x + tmp);
+ sin_lamino = sincos (mad((float) idz, lamino_region.y, lamino_region.x), &cos_lamino);
+ tmp = mad(z_region.x + 0.5f + ((int) y_center) - y_center, sin_lamino, y_center);
+ tmp_x = voxel.x * cos_lamino;
+ tmp_y = -voxel.y * cos_lamino;
+
+{}
diff --git a/src/kernels/templates/roll_template.in b/src/kernels/templates/roll_template.in
new file mode 100644
index 0000000..67addea
--- /dev/null
+++ b/src/kernels/templates/roll_template.in
@@ -0,0 +1,14 @@
+ pixel.x = mad(voxel.x, cosines{1}, mad(voxel.y, sines{1}, x_center.x));
+ pixel.y = mad(tmp_x, sines{1}, mad(tmp_y, cosines{1}, tmp));
+ rotate ();
+ result {2}= read_imagef (projection_{0}, sampler, pixel).x;
+%nl
+ tmp = 0.5f + ((int) x_center.x) - x_center.x;
+ voxel.x = mad((float) idx, x_region.y, x_region.x + tmp);
+ voxel.y = mad((float) idy, y_region.y, y_region.x + tmp);
+ sin_roll = sincos (mad((float) idz, -roll_region.y, -roll_region.x), &cos_roll);
+ tmp = mad(z_region.x + 0.5f + ((int) y_center) - y_center, sin_lamino, y_center);
+ tmp_x = voxel.x * cos_lamino;
+ tmp_y = -voxel.y * cos_lamino;
+
+{}
diff --git a/src/kernels/templates/z_template.in b/src/kernels/templates/z_template.in
new file mode 100644
index 0000000..a0d46c0
--- /dev/null
+++ b/src/kernels/templates/z_template.in
@@ -0,0 +1,19 @@
+ pixel.x = mad(voxel.x, cosines{1}, mad(voxel.y, sines{1}, x_center.x));
+ pixel.y = mad(tmp_x, sines{1}, mad(tmp_y, cosines{1}, tmp));
+ rotate ();
+ result {2}= read_imagef (projection_{0}, sampler, pixel).x;
+%nl
+ /* Pixel interval is (0, 1), so it is evaluated at 0.5. The volume is */
+ /* shifted so that the pixel is rotated with the +0.5 offset if the axis */
+ /* falls between two pixels and with 0 offset if the axis is in the middle */
+ /* of a pixel */
+ /* Simple arithmetic faster than modf */
+ tmp = 0.5f + ((int) x_center.x) - x_center.x;
+ voxel.x = mad((float) idx, x_region.y, x_region.x + tmp);
+ voxel.y = mad((float) idy, y_region.y, y_region.x + tmp);
+ voxel.z = mad((float) idz, z_region.y, z_region.x + 0.5f + ((int) y_center) - y_center);
+ tmp = mad(voxel.z, sin_lamino, y_center);
+ tmp_x = voxel.x * cos_lamino;
+ tmp_y = -voxel.y * cos_lamino;
+
+{}
diff --git a/src/kernels/tools/make_burst_kernels.py b/src/kernels/tools/make_burst_kernels.py
new file mode 100644
index 0000000..d403c93
--- /dev/null
+++ b/src/kernels/tools/make_burst_kernels.py
@@ -0,0 +1,63 @@
+"""Generate burst laminographic backprojection OpenCL kernels."""
+import argparse
+import os
+
+
+IDX_TO_VEC_ELEM = dict(zip(range(10), range(10)))
+IDX_TO_VEC_ELEM[10] = 'a'
+IDX_TO_VEC_ELEM[11] = 'b'
+IDX_TO_VEC_ELEM[12] = 'c'
+IDX_TO_VEC_ELEM[13] = 'd'
+IDX_TO_VEC_ELEM[14] = 'e'
+IDX_TO_VEC_ELEM[15] = 'f'
+
+
+def fill_compute_template(tmpl, num_items, index):
+ """Fill the template doing the pixel computation and texture fetch."""
+ operation = '+' if index else ''
+ access = '.s{}'.format(IDX_TO_VEC_ELEM[index]) if num_items > 1 else ''
+ return tmpl.format(index, access, operation)
+
+
+def fill_kernel_template(input_tmpl, compute_tmpl, kernel_outer, kernel_inner, num_items):
+ """Construct the whole kernel."""
+ vector_length = num_items if num_items > 1 else ''
+ computes = '\n'.join([fill_compute_template(compute_tmpl, num_items, i)
+ for i in range(num_items)])
+ inputs = '\n'.join([input_tmpl.format(i) for i in range(num_items)])
+ kernel_inner = kernel_inner.format(computes)
+
+ return kernel_outer.format(num_items, inputs, vector_length, kernel_inner)
+
+
+def parse_args():
+ """Parse command line arguments."""
+ parser = argparse.ArgumentParser()
+ parser.add_argument('filename', type=str, help='File name with the kernel template')
+ parser.add_argument('burst', type=int, nargs='+',
+ help='Number of projections processed by one kernel invocation')
+
+ return parser.parse_args()
+
+
+def main():
+ """execute program."""
+ args = parse_args()
+ allowed_bursts = [2 ** i for i in range(5)]
+ in_tmpl = "read_only image2d_t projection_{},"
+ common_filename = os.path.join(os.path.dirname(args.filename), 'common.in')
+ defs_filename = os.path.join(os.path.dirname(args.filename), 'definitions.in')
+ defs = open(defs_filename, 'r').read()
+ kernel_outer = open(common_filename, 'r').read()
+ comp_tmpl, kernel_inner = open(args.filename, 'r').read().split('\n%nl\n')
+ kernels = defs + '\n'
+ for burst in args.burst:
+ if burst not in allowed_bursts:
+ raise ValueError('Burst mode `{}` must be one of `{}`'.format(burst, allowed_bursts))
+ kernels += fill_kernel_template(in_tmpl, comp_tmpl, kernel_outer, kernel_inner, burst)
+
+ print kernels
+
+
+if __name__ == '__main__':
+ main()
diff --git a/src/lamino-roi.c b/src/lamino-roi.c
new file mode 100644
index 0000000..c7875a9
--- /dev/null
+++ b/src/lamino-roi.c
@@ -0,0 +1,153 @@
+#include <math.h>
+#include <glib.h>
+#include "lamino-roi.h"
+
+
+static inline void
+swap (gint *first, gint *second) {
+ gint tmp;
+
+ tmp = *first;
+ *first = *second;
+ *second = tmp;
+}
+
+/**
+ * Determine the left and right column to read from a projection at a given
+ * tomographic angle.
+ */
+void
+determine_x_extrema (gfloat extrema[2], GValueArray *x_extrema, GValueArray *y_extrema,
+ gfloat tomo_angle, gfloat x_center)
+{
+ gfloat sin_tomo, cos_tomo;
+ gint x_min, x_max, y_min, y_max;
+
+ sin_tomo = sin (tomo_angle);
+ cos_tomo = cos (tomo_angle);
+ x_min = EXTRACT_INT (x_extrema, 0);
+ /* The interval is right-opened when OpenCL indices for both x and y are generated, */
+ /* so the last index doesn't count */
+ x_max = EXTRACT_INT (x_extrema, 1) - 1;
+ y_min = EXTRACT_INT (y_extrema, 0);
+ y_max = EXTRACT_INT (y_extrema, 1) - 1;
+
+ if (sin_tomo < 0) {
+ swap (&y_min, &y_max);
+ }
+ if (cos_tomo < 0) {
+ swap (&x_min, &x_max);
+ }
+
+ /* -1 to make sure the interpolation doesn't reach to uninitialized values*/
+ extrema[0] = cos_tomo * x_min + sin_tomo * y_min + x_center - 1;
+ /* +1 because extrema[1] will be accessed by interpolation
+ * but the region in copying is right-open */
+ extrema[1] = cos_tomo * x_max + sin_tomo * y_max + x_center + 1;
+}
+
+/**
+ * Determine the top and bottom row to read from a projection at given
+ * tomographic and laminographic angles.
+ */
+void
+determine_y_extrema (gfloat extrema[2], GValueArray *x_extrema, GValueArray *y_extrema,
+ gfloat z_extrema[2], gfloat tomo_angle, gfloat lamino_angle,
+ gfloat y_center)
+{
+ gfloat sin_tomo, cos_tomo, sin_lamino, cos_lamino;
+ gint x_min, x_max, y_min, y_max;
+
+ sin_tomo = sin (tomo_angle);
+ cos_tomo = cos (tomo_angle);
+ sin_lamino = sin (lamino_angle);
+ cos_lamino = cos (lamino_angle);
+ x_min = EXTRACT_INT (x_extrema, 0);
+ x_max = EXTRACT_INT (x_extrema, 1) - 1;
+ y_min = EXTRACT_INT (y_extrema, 0);
+ y_max = EXTRACT_INT (y_extrema, 1) - 1;
+
+ if (sin_tomo < 0) {
+ swap (&x_min, &x_max);
+ }
+ if (cos_tomo > 0) {
+ swap (&y_min, &y_max);
+ }
+
+ extrema[0] = sin_tomo * x_min - cos_tomo * y_min;
+ extrema[1] = sin_tomo * x_max - cos_tomo * y_max;
+ extrema[0] = extrema[0] * cos_lamino + z_extrema[0] * sin_lamino + y_center - 1;
+ extrema[1] = extrema[1] * cos_lamino + z_extrema[1] * sin_lamino + y_center + 1;
+}
+
+/**
+ * clip:
+ * @result: resulting clipped extrema
+ * @extrema: (min, max)
+ * @maximum: projection width or height
+ *
+ * Clip extrema to an allowed interval [0, projection width/height)
+ */
+void
+clip (gint result[2], gfloat extrema[2], gint maximum)
+{
+ result[0] = (gint) floorf (extrema[0]);
+ result[1] = (gint) ceilf (extrema[1]);
+
+ if (result[0] < 0) {
+ result[0] = 0;
+ }
+ if (result[0] > maximum) {
+ result[0] = maximum;
+ }
+ if (result[1] < 0) {
+ result[1] = 0;
+ }
+ if (result[1] > maximum) {
+ result[1] = maximum;
+ }
+
+ if (result[0] == result[1]) {
+ if (result[1] < maximum) {
+ result[1]++;
+ } else if (result[0] > 0) {
+ result[0]--;
+ } else {
+ g_warning ("Cannot extend");
+ }
+ } else if (result[0] > result[1]) {
+ g_warning ("Invalid extrema: minimum larger than maximum");
+ }
+}
+
+/**
+ * Determine the left and right column to read from a projection at a given
+ * tomographic angle. The result is bound to [0, projection width)
+ */
+void
+determine_x_region (gint result[2], GValueArray *x_extrema, GValueArray *y_extrema, gfloat tomo_angle,
+ gfloat x_center, gint width)
+{
+ gfloat extrema[2];
+
+ determine_x_extrema (extrema, x_extrema, y_extrema, tomo_angle, x_center);
+
+ clip (result, extrema, width);
+}
+
+/**
+ * Determine the top and bottom column to read from a projection at given
+ * tomographic and laminographic angles. The result is bound to
+ * [0, projection height).
+ */
+void
+determine_y_region (gint result[2], GValueArray *x_extrema, GValueArray *y_extrema, gfloat z_extrema[2],
+ gfloat tomo_angle, gfloat lamino_angle, gfloat y_center, gint height)
+{
+ gfloat extrema[2];
+
+ determine_y_extrema (extrema, x_extrema, y_extrema, z_extrema, tomo_angle,
+ lamino_angle, y_center);
+ clip (result, extrema, height);
+}
+
diff --git a/src/lamino-roi.h b/src/lamino-roi.h
new file mode 100644
index 0000000..ff6876d
--- /dev/null
+++ b/src/lamino-roi.h
@@ -0,0 +1,23 @@
+#ifndef LAMINO_ROI_H
+#define LAMINO_ROI_H
+
+#define EXTRACT_INT(region, index) g_value_get_int (g_value_array_get_nth ((region), (index)))
+
+#include <glib-object.h>
+
+G_BEGIN_DECLS
+
+void clip (gint result[2], gfloat extrema[2], gint maximum);
+void determine_x_extrema (gfloat extrema[2], GValueArray *x_extrema, GValueArray *y_extrema,
+ gfloat tomo_angle, gfloat x_center);
+void determine_y_extrema (gfloat extrema[2], GValueArray *x_extrema, GValueArray *y_extrema,
+ gfloat z_extrema[2], gfloat tomo_angle, gfloat lamino_angle,
+ gfloat y_center);
+void determine_x_region (gint result[2], GValueArray *x_extrema, GValueArray *y_extrema, gfloat tomo_angle,
+ gfloat x_center, gint width);
+void determine_y_region (gint result[2], GValueArray *x_extrema, GValueArray *y_extrema, gfloat z_extrema[2],
+ gfloat tomo_angle, gfloat lamino_angle, gfloat y_center, gint height);
+
+G_END_DECLS
+
+#endif
diff --git a/src/ufo-lamino-backproject-task.c b/src/ufo-lamino-backproject-task.c
new file mode 100644
index 0000000..bf74271
--- /dev/null
+++ b/src/ufo-lamino-backproject-task.c
@@ -0,0 +1,828 @@
+/*
+ * Copyright (C) 2011-2015 Karlsruhe Institute of Technology
+ *
+ * This file is part of Ufo.
+ *
+ * This library is free software: you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public
+ * License as published by the Free Software Foundation, either
+ * version 3 of the License, or (at your option) any later version.
+ *
+ * This library is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this library. If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#include <stdio.h>
+#include <math.h>
+#include <glib.h>
+#include <glib/gprintf.h>
+
+#ifdef __APPLE__
+#include <OpenCL/cl.h>
+#else
+#include <CL/cl.h>
+#endif
+
+#include "ufo-lamino-backproject-task.h"
+#include "lamino-roi.h"
+
+/* Copy only neccessary projection region */
+/* TODO: make this a parameter? */
+/* Wait with enabling this until sync issues in ufo-core have been solved */
+#define COPY_PROJECTION_REGION 0
+#define EXTRACT_FLOAT(region, index) g_value_get_float (g_value_array_get_nth ((region), (index)))
+#define REGION_SIZE(region) ((EXTRACT_INT ((region), 2)) == 0) ? 0 : \
+ ((EXTRACT_INT ((region), 1) - EXTRACT_INT ((region), 0) - 1) /\
+ EXTRACT_INT ((region), 2) + 1)
+#define PAD_TO_DIVIDE(dividend, divisor) ((dividend) + (divisor) - (dividend) % (divisor))
+
+
+typedef enum {
+ PARAM_Z,
+ PARAM_CENTER,
+ PARAM_LAMINO,
+ PARAM_ROLL
+} Param;
+
+struct _UfoLaminoBackprojectTaskPrivate {
+ /* private */
+ gboolean generated;
+ guint count;
+ /* sine and cosine table size based on BURST */
+ gsize table_size;
+
+ /* OpenCL */
+ cl_context context;
+ cl_kernel vector_kernel;
+ cl_kernel scalar_kernel;
+ cl_sampler sampler;
+ /* Buffered images for invoking backprojection on BURST projections at once.
+ * We potentially don't need to copy the last image and can use the one from
+ * framework directly but it seems to have no performance effects. */
+ cl_mem images[BURST];
+
+ /* properties */
+ GValueArray *x_region;
+ GValueArray *y_region;
+ GValueArray *region;
+ GValueArray *center;
+ GValueArray *projection_offset;
+ float sines[BURST], cosines[BURST];
+ guint num_projections;
+ gfloat overall_angle;
+ gfloat tomo_angle;
+ gfloat lamino_angle;
+ gfloat z;
+ gfloat roll_angle;
+ Param parameter;
+};
+
+static void ufo_task_interface_init (UfoTaskIface *iface);
+
+G_DEFINE_TYPE_WITH_CODE (UfoLaminoBackprojectTask, ufo_lamino_backproject_task, UFO_TYPE_TASK_NODE,
+ G_IMPLEMENT_INTERFACE (UFO_TYPE_TASK,
+ ufo_task_interface_init))
+
+#define UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE(obj) (G_TYPE_INSTANCE_GET_PRIVATE((obj), UFO_TYPE_LAMINO_BACKPROJECT_TASK, UfoLaminoBackprojectTaskPrivate))
+
+enum {
+ PROP_0,
+ PROP_X_REGION,
+ PROP_Y_REGION,
+ PROP_Z,
+ PROP_REGION,
+ PROP_PROJECTION_OFFSET,
+ PROP_CENTER,
+ PROP_NUM_PROJECTIONS,
+ PROP_OVERALL_ANGLE,
+ PROP_TOMO_ANGLE,
+ PROP_LAMINO_ANGLE,
+ PROP_PARAMETER,
+ PROP_ROLL_ANGLE,
+ N_PROPERTIES
+};
+
+static GParamSpec *properties[N_PROPERTIES] = { NULL, };
+
+static void
+set_region (GValueArray *src, GValueArray **dst)
+{
+ if (EXTRACT_INT (src, 0) > EXTRACT_INT (src, 1)) {
+ g_warning ("Invalid region [\"from\", \"to\", \"step\"]: [%d, %d, %d], "\
+ "\"from\" has to be less than or equal to \"to\"",
+ EXTRACT_INT (src, 0), EXTRACT_INT (src, 1), EXTRACT_INT (src, 2));
+ }
+ else {
+ g_value_array_free (*dst);
+ *dst = g_value_array_copy (src);
+ }
+}
+
+static void
+copy_to_image (UfoBuffer *input,
+ cl_mem output_image,
+ cl_command_queue cmd_queue,
+ size_t origin[3],
+ size_t region[3],
+ gint in_width)
+{
+ cl_mem input_data;
+ cl_event event;
+
+ input_data = ufo_buffer_get_device_image (input, cmd_queue);
+ UFO_RESOURCES_CHECK_CLERR (clEnqueueCopyImage (cmd_queue,
+ input_data,
+ output_image,
+ origin,
+ origin,
+ region,
+ 0,
+ NULL,
+ &event));
+
+ UFO_RESOURCES_CHECK_CLERR (clWaitForEvents (1, &event));
+ UFO_RESOURCES_CHECK_CLERR (clReleaseEvent (event));
+}
+
+UfoNode *
+ufo_lamino_backproject_task_new (void)
+{
+ return UFO_NODE (g_object_new (UFO_TYPE_LAMINO_BACKPROJECT_TASK, NULL));
+}
+
+static void
+ufo_lamino_backproject_task_setup (UfoTask *task,
+ UfoResources *resources,
+ GError **error)
+{
+ UfoLaminoBackprojectTaskPrivate *priv;
+ cl_int cl_error;
+ gint i;
+ gchar *vector_kernel_name, *kernel_filename;
+
+ vector_kernel_name = g_strdup_printf ("backproject_burst_%d", BURST);
+ if (!vector_kernel_name) {
+ g_warning ("Error making burst kernel name");
+ }
+
+ priv = UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE (task);
+ priv->context = ufo_resources_get_context (resources);
+ switch (priv->parameter) {
+ case PARAM_Z:
+ kernel_filename = g_strdup ("z_kernel.cl");
+ break;
+ case PARAM_CENTER:
+ kernel_filename = g_strdup ("center_kernel.cl");
+ break;
+ case PARAM_LAMINO:
+ kernel_filename = g_strdup ("lamino_kernel.cl");
+ break;
+ case PARAM_ROLL:
+ kernel_filename = g_strdup ("roll_kernel.cl");
+ break;
+ default:
+ g_warning ("Unkown varying parameter");
+ break;
+ }
+ priv->vector_kernel = ufo_resources_get_kernel (resources, kernel_filename,
+ vector_kernel_name, error);
+ priv->scalar_kernel = ufo_resources_get_kernel (resources, kernel_filename,
+ "backproject_burst_1", error);
+ priv->sampler = clCreateSampler (priv->context,
+ (cl_bool) FALSE,
+ CL_ADDRESS_CLAMP,
+ CL_FILTER_LINEAR,
+ &cl_error);
+
+ UFO_RESOURCES_CHECK_CLERR (clRetainContext (priv->context));
+ UFO_RESOURCES_CHECK_CLERR (cl_error);
+ if (priv->vector_kernel) {
+ UFO_RESOURCES_CHECK_CLERR (clRetainKernel (priv->vector_kernel));
+ }
+ if (priv->scalar_kernel) {
+ UFO_RESOURCES_CHECK_CLERR (clRetainKernel (priv->scalar_kernel));
+ }
+
+ for (i = 0; i < BURST; i++) {
+ priv->images[i] = NULL;
+ }
+
+ switch (BURST) {
+ case 1: priv->table_size = sizeof (cl_float); break;
+ case 2: priv->table_size = sizeof (cl_float2); break;
+ case 4: priv->table_size = sizeof (cl_float4); break;
+ case 8: priv->table_size = sizeof (cl_float8); break;
+ case 16: priv->table_size = sizeof (cl_float16); break;
+ default: g_warning ("Unsupported vector size"); break;
+ }
+
+ g_free (vector_kernel_name);
+ g_free (kernel_filename);
+}
+
+static void
+ufo_lamino_backproject_task_get_requisition (UfoTask *task,
+ UfoBuffer **inputs,
+ UfoRequisition *requisition)
+{
+ UfoLaminoBackprojectTaskPrivate *priv;
+ gfloat start, stop, step;
+
+ priv = UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE (task);
+ start = EXTRACT_FLOAT (priv->region, 0);
+ stop = EXTRACT_FLOAT (priv->region, 1);
+ step = EXTRACT_FLOAT (priv->region, 2);
+
+ if (!priv->num_projections) {
+ g_warning ("Number of projections has not been set");
+ }
+ if (step == 0.0f) {
+ g_warning ("Step in region is 0");
+ }
+
+ requisition->n_dims = 3;
+ requisition->dims[0] = REGION_SIZE (priv->x_region);
+ requisition->dims[1] = REGION_SIZE (priv->y_region);
+ requisition->dims[2] = (gint) ceil ((stop - start) / step);
+}
+
+static guint
+ufo_lamino_backproject_task_get_num_inputs (UfoTask *task)
+{
+ return 1;
+}
+
+static guint
+ufo_lamino_backproject_task_get_num_dimensions (UfoTask *task,
+ guint input)
+{
+ g_return_val_if_fail (input == 0, 0);
+
+ return 3;
+}
+
+static UfoTaskMode
+ufo_lamino_backproject_task_get_mode (UfoTask *task)
+{
+ return UFO_TASK_MODE_REDUCTOR | UFO_TASK_MODE_GPU;
+}
+
+static gboolean
+ufo_lamino_backproject_task_process (UfoTask *task,
+ UfoBuffer **inputs,
+ UfoBuffer *output,
+ UfoRequisition *requisition)
+{
+ UfoLaminoBackprojectTaskPrivate *priv;
+ UfoRequisition in_req;
+ UfoGpuNode *node;
+ UfoProfiler *profiler;
+ gfloat tomo_angle, *sines, *cosines;
+ gint i, index;
+ gint cumulate;
+ gsize table_size;
+ gboolean scalar;
+ /* regions stripped off the "to" value */
+ gfloat x_region[2], y_region[2], z_region[2], x_center[2], z_ends[2], lamino_angles[2], roll_angles[2],
+ y_center, sin_lamino, cos_lamino, norm_factor, sin_roll, cos_roll;
+ gint x_copy_region[2], y_copy_region[2];
+ cl_kernel kernel;
+ cl_command_queue cmd_queue;
+ cl_mem out_mem;
+ cl_int cl_error;
+ /* image creation and copying */
+ cl_image_format image_fmt;
+ size_t origin[3];
+ size_t region[3];
+ /* keep the warp size satisfied but make sure the local grid is localized
+ * around a point in 3D for efficient caching */
+ const gint real_size[4] = {requisition->dims[0], requisition->dims[1], requisition->dims[2], 0};
+ const gsize local_work_size[] = {16, 8, 8};
+ gsize global_work_size[3];
+
+ global_work_size[0] = requisition->dims[0] % local_work_size[0] ?
+ PAD_TO_DIVIDE (requisition->dims[0], local_work_size[0]) :
+ requisition->dims[0];
+ global_work_size[1] = requisition->dims[1] % local_work_size[1] ?
+ PAD_TO_DIVIDE (requisition->dims[1], local_work_size[1]) :
+ requisition->dims[1];
+ global_work_size[2] = requisition->dims[2] % local_work_size[2] ?
+ PAD_TO_DIVIDE (requisition->dims[2], local_work_size[2]) :
+ requisition->dims[2];
+
+ priv = UFO_LAMINO_BACKPROJECT_TASK (task)->priv;
+ node = UFO_GPU_NODE (ufo_task_node_get_proc_node (UFO_TASK_NODE (task)));
+ cmd_queue = ufo_gpu_node_get_cmd_queue (node);
+ out_mem = ufo_buffer_get_device_array (output, cmd_queue);
+ ufo_buffer_get_requisition (inputs[0], &in_req);
+
+ index = priv->count % BURST;
+ tomo_angle = priv->tomo_angle > -G_MAXFLOAT ? priv->tomo_angle :
+ priv->overall_angle * priv->count / priv->num_projections;
+ norm_factor = fabs (priv->overall_angle) / priv->num_projections;
+ priv->sines[index] = sin (tomo_angle);
+ priv->cosines[index] = cos (tomo_angle);
+ x_region[0] = (gfloat) EXTRACT_INT (priv->x_region, 0);
+ x_region[1] = (gfloat) EXTRACT_INT (priv->x_region, 2);
+ y_region[0] = (gfloat) EXTRACT_INT (priv->y_region, 0);
+ y_region[1] = (gfloat) EXTRACT_INT (priv->y_region, 2);
+ if (priv->parameter == PARAM_Z) {
+ z_ends[0] = z_region[0] = EXTRACT_FLOAT (priv->region, 0);
+ z_region[1] = EXTRACT_FLOAT (priv->region, 2);
+ z_ends[1] = EXTRACT_FLOAT (priv->region, 1);
+ } else {
+ z_ends[0] = z_region[0] = priv->z;
+ z_ends[1] = priv->z + 1.0f;
+ }
+ if (priv->parameter == PARAM_CENTER) {
+ x_center[0] = EXTRACT_FLOAT (priv->region, 0) - EXTRACT_INT (priv->projection_offset, 0);
+ x_center[1] = EXTRACT_FLOAT (priv->region, 2);
+ } else {
+ x_center[0] = x_center[1] = EXTRACT_FLOAT (priv->center, 0) - EXTRACT_INT (priv->projection_offset, 0);
+ }
+ if (priv->parameter == PARAM_LAMINO) {
+ lamino_angles[0] = EXTRACT_FLOAT (priv->region, 0);
+ lamino_angles[1] = EXTRACT_FLOAT (priv->region, 2);
+ } else {
+ lamino_angles[0] = lamino_angles[1] = priv->lamino_angle;
+ }
+ if (priv->parameter == PARAM_ROLL) {
+ roll_angles[0] = EXTRACT_FLOAT (priv->region, 0);
+ roll_angles[1] = EXTRACT_FLOAT (priv->region, 2);
+ } else {
+ roll_angles[0] = roll_angles[1] = priv->roll_angle;
+ }
+
+ y_center = EXTRACT_FLOAT (priv->center, 1) - EXTRACT_INT (priv->projection_offset, 1);
+ sin_lamino = sinf (priv->lamino_angle);
+ cos_lamino = cosf (priv->lamino_angle);
+ /* Minus the value because we are rotating back */
+ sin_roll = sinf (-priv->roll_angle);
+ cos_roll = cosf (-priv->roll_angle);
+ scalar = priv->count >= priv->num_projections / BURST * BURST ? 1 : 0;
+
+ /* If COPY_PROJECTION_REGION is True we copy only the part necessary */
+ /* for a given tomographic and laminographic angle */
+ /* TODO: Extend the region determination to be able to handle PARAM_LAMINO */
+ if (COPY_PROJECTION_REGION && priv->parameter != PARAM_LAMINO) {
+ determine_x_region (x_copy_region, priv->x_region, priv->y_region, tomo_angle,
+ EXTRACT_FLOAT (priv->center, 0), in_req.dims[0]);
+ determine_y_region (y_copy_region, priv->x_region, priv->y_region, z_ends,
+ tomo_angle, priv->lamino_angle, EXTRACT_FLOAT (priv->center, 1),
+ in_req.dims[1]);
+ origin[0] = x_copy_region[0];
+ origin[1] = y_copy_region[0];
+ origin[2] = 0;
+ region[0] = x_copy_region[1] - x_copy_region[0];
+ region[1] = y_copy_region[1] - y_copy_region[0];
+ } else {
+ origin[0] = origin[1] = origin[2] = 0;
+ region[0] = in_req.dims[0];
+ region[1] = in_req.dims[1];
+ }
+ region[2] = 1;
+
+ if (priv->images[index] == NULL) {
+ /* TODO: dangerous, don't rely on the ufo-buffer */
+ image_fmt.image_channel_order = CL_INTENSITY;
+ image_fmt.image_channel_data_type = CL_FLOAT;
+ /* TODO: what with the "other" API? */
+ priv->images[index] = clCreateImage2D (priv->context,
+ CL_MEM_READ_ONLY,
+ &image_fmt,
+ in_req.dims[0],
+ in_req.dims[1],
+ 0,
+ NULL,
+ &cl_error);
+ UFO_RESOURCES_CHECK_CLERR (cl_error);
+ }
+
+ copy_to_image (inputs[0], priv->images[index], cmd_queue, origin, region, in_req.dims[0]);
+
+ if (scalar) {
+ kernel = priv->scalar_kernel;
+ cumulate = priv->count;
+ table_size = sizeof (cl_float);
+ sines = &priv->sines[index];
+ cosines = &priv->cosines[index];
+ i = 1;
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, 0, sizeof (cl_mem), &priv->images[index]));
+ } else {
+ kernel = priv->vector_kernel;
+ cumulate = priv->count + 1 == BURST ? 0 : 1;
+ table_size = priv->table_size;
+ sines = priv->sines;
+ cosines = priv->cosines;
+ i = BURST;
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, index, sizeof (cl_mem), &priv->images[index]));
+ }
+
+ if (scalar || index == BURST - 1) {
+ /* Execute the kernel after BURST images have arrived, i.e. we use more
+ * projections at one invocation, so the number of read/writes to the
+ * result is reduced by a factor of BURST. If there are not enough
+ * projecttions left, execute the scalar kernel */
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_mem), &out_mem));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_sampler), &priv->sampler));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_int3), real_size));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float2), x_center));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float), (cl_float *) &y_center));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float2), x_region));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float2), y_region));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float2), z_region));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float2), lamino_angles));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float2), roll_angles));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float), &sin_lamino));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float), &cos_lamino));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, table_size, sines));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, table_size, cosines));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float), &norm_factor));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float), &sin_roll));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i++, sizeof (cl_float), &cos_roll));
+ UFO_RESOURCES_CHECK_CLERR (clSetKernelArg (kernel, i, sizeof (cl_int), (cl_int *) &cumulate));
+
+ profiler = ufo_task_node_get_profiler (UFO_TASK_NODE (task));
+ ufo_profiler_call (profiler, cmd_queue, kernel, 3, global_work_size, local_work_size);
+ }
+
+ priv->count++;
+
+ return TRUE;
+}
+
+static gboolean
+ufo_lamino_backproject_task_generate (UfoTask *task,
+ UfoBuffer *output,
+ UfoRequisition *requisition)
+{
+ UfoLaminoBackprojectTaskPrivate *priv;
+
+ priv = UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE (task);
+
+ if (priv->generated) {
+ return FALSE;
+ }
+
+ priv->generated = TRUE;
+
+ return TRUE;
+}
+
+static void
+ufo_lamino_backproject_task_set_property (GObject *object,
+ guint property_id,
+ const GValue *value,
+ GParamSpec *pspec)
+{
+ UfoLaminoBackprojectTaskPrivate *priv = UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE (object);
+ GValueArray *array;
+
+ switch (property_id) {
+ case PROP_X_REGION:
+ array = (GValueArray *) g_value_get_boxed (value);
+ set_region (array, &priv->x_region);
+ break;
+ case PROP_Y_REGION:
+ array = (GValueArray *) g_value_get_boxed (value);
+ set_region (array, &priv->y_region);
+ break;
+ case PROP_Z:
+ priv->z = g_value_get_float (value);
+ break;
+ case PROP_REGION:
+ array = (GValueArray *) g_value_get_boxed (value);
+ g_value_array_free (priv->region);
+ priv->region = g_value_array_copy (array);
+ break;
+ case PROP_PROJECTION_OFFSET:
+ array = (GValueArray *) g_value_get_boxed (value);
+ g_value_array_free (priv->projection_offset);
+ priv->projection_offset = g_value_array_copy (array);
+ break;
+ case PROP_CENTER:
+ array = (GValueArray *) g_value_get_boxed (value);
+ g_value_array_free (priv->center);
+ priv->center = g_value_array_copy (array);
+ break;
+ case PROP_NUM_PROJECTIONS:
+ priv->num_projections = g_value_get_uint (value);
+ break;
+ case PROP_OVERALL_ANGLE:
+ priv->overall_angle = g_value_get_float (value);
+ break;
+ case PROP_TOMO_ANGLE:
+ priv->tomo_angle = g_value_get_float (value);
+ break;
+ case PROP_LAMINO_ANGLE:
+ priv->lamino_angle = g_value_get_float (value);
+ break;
+ case PROP_ROLL_ANGLE:
+ priv->roll_angle = g_value_get_float (value);
+ break;
+ case PROP_PARAMETER:
+ if (!g_strcmp0 (g_value_get_string (value), "z")) {
+ priv->parameter = PARAM_Z;
+ } else if (!g_strcmp0 (g_value_get_string (value), "x-center")) {
+ priv->parameter = PARAM_CENTER;
+ } else if (!g_strcmp0 (g_value_get_string (value), "lamino-angle")) {
+ priv->parameter = PARAM_LAMINO;
+ } else if (!g_strcmp0 (g_value_get_string (value), "roll-angle")) {
+ priv->parameter = PARAM_ROLL;
+ }
+ break;
+ default:
+ G_OBJECT_WARN_INVALID_PROPERTY_ID (object, property_id, pspec);
+ break;
+ }
+}
+
+static void
+ufo_lamino_backproject_task_get_property (GObject *object,
+ guint property_id,
+ GValue *value,
+ GParamSpec *pspec)
+{
+ UfoLaminoBackprojectTaskPrivate *priv = UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE (object);
+
+ switch (property_id) {
+ case PROP_X_REGION:
+ g_value_set_boxed (value, priv->x_region);
+ break;
+ case PROP_Y_REGION:
+ g_value_set_boxed (value, priv->y_region);
+ break;
+ case PROP_Z:
+ g_value_set_float (value, priv->z);
+ break;
+ case PROP_REGION:
+ g_value_set_boxed (value, priv->region);
+ break;
+ case PROP_PROJECTION_OFFSET:
+ g_value_set_boxed (value, priv->projection_offset);
+ break;
+ case PROP_CENTER:
+ g_value_set_boxed (value, priv->center);
+ break;
+ case PROP_NUM_PROJECTIONS:
+ g_value_set_uint (value, priv->num_projections);
+ break;
+ case PROP_OVERALL_ANGLE:
+ g_value_set_float (value, priv->overall_angle);
+ break;
+ case PROP_TOMO_ANGLE:
+ g_value_set_float (value, priv->tomo_angle);
+ break;
+ case PROP_LAMINO_ANGLE:
+ g_value_set_float (value, priv->lamino_angle);
+ break;
+ case PROP_ROLL_ANGLE:
+ g_value_set_float (value, priv->roll_angle);
+ break;
+ case PROP_PARAMETER:
+ switch (priv->parameter) {
+ case PARAM_Z:
+ g_value_set_string (value, "z");
+ break;
+ case PARAM_CENTER:
+ g_value_set_string (value, "x-center");
+ break;
+ case PARAM_LAMINO:
+ g_value_set_string (value, "lamino-angle");
+ break;
+ case PARAM_ROLL:
+ g_value_set_string (value, "roll-angle");
+ break;
+ }
+ break;
+ default:
+ G_OBJECT_WARN_INVALID_PROPERTY_ID (object, property_id, pspec);
+ break;
+ }
+}
+
+static void
+ufo_lamino_backproject_task_finalize (GObject *object)
+{
+ UfoLaminoBackprojectTaskPrivate *priv;
+ gint i;
+
+ priv = UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE (object);
+ g_value_array_free (priv->x_region);
+ g_value_array_free (priv->y_region);
+ g_value_array_free (priv->region);
+ g_value_array_free (priv->projection_offset);
+ g_value_array_free (priv->center);
+
+ if (priv->vector_kernel) {
+ UFO_RESOURCES_CHECK_CLERR (clReleaseKernel (priv->vector_kernel));
+ priv->vector_kernel = NULL;
+ }
+ if (priv->scalar_kernel) {
+ UFO_RESOURCES_CHECK_CLERR (clReleaseKernel (priv->scalar_kernel));
+ priv->scalar_kernel = NULL;
+ }
+ if (priv->context) {
+ UFO_RESOURCES_CHECK_CLERR (clReleaseContext (priv->context));
+ priv->context = NULL;
+ }
+ if (priv->sampler) {
+ UFO_RESOURCES_CHECK_CLERR (clReleaseSampler (priv->sampler));
+ priv->sampler = NULL;
+ }
+
+ for (i = 0; i < BURST; i++) {
+ if (priv->images[i] != NULL) {
+ UFO_RESOURCES_CHECK_CLERR (clReleaseMemObject (priv->images[i]));
+ priv->images[i] = NULL;
+ }
+ }
+
+ G_OBJECT_CLASS (ufo_lamino_backproject_task_parent_class)->finalize (object);
+}
+
+static void
+ufo_task_interface_init (UfoTaskIface *iface)
+{
+ iface->setup = ufo_lamino_backproject_task_setup;
+ iface->get_num_inputs = ufo_lamino_backproject_task_get_num_inputs;
+ iface->get_num_dimensions = ufo_lamino_backproject_task_get_num_dimensions;
+ iface->get_mode = ufo_lamino_backproject_task_get_mode;
+ iface->get_requisition = ufo_lamino_backproject_task_get_requisition;
+ iface->process = ufo_lamino_backproject_task_process;
+ iface->generate = ufo_lamino_backproject_task_generate;
+}
+
+static void
+ufo_lamino_backproject_task_class_init (UfoLaminoBackprojectTaskClass *klass)
+{
+ GObjectClass *oclass = G_OBJECT_CLASS (klass);
+
+ oclass->set_property = ufo_lamino_backproject_task_set_property;
+ oclass->get_property = ufo_lamino_backproject_task_get_property;
+ oclass->finalize = ufo_lamino_backproject_task_finalize;
+
+ GParamSpec *region_vals = g_param_spec_int ("region_values",
+ "Region values",
+ "Elements in regions",
+ G_MININT,
+ G_MAXINT,
+ (gint) 0,
+ G_PARAM_READWRITE);
+
+ GParamSpec *float_region_vals = g_param_spec_float ("float_region_values",
+ "Float Region values",
+ "Elements in float regions",
+ -G_MAXFLOAT,
+ G_MAXFLOAT,
+ 0.0f,
+ G_PARAM_READWRITE);
+
+ properties[PROP_X_REGION] =
+ g_param_spec_value_array ("x-region",
+ "X region for reconstruction as (from, to, step)",
+ "X region for reconstruction as (from, to, step)",
+ region_vals,
+ G_PARAM_READWRITE);
+
+ properties[PROP_Y_REGION] =
+ g_param_spec_value_array ("y-region",
+ "Y region for reconstruction as (from, to, step)",
+ "Y region for reconstruction as (from, to, step)",
+ region_vals,
+ G_PARAM_READWRITE);
+
+ properties[PROP_Z] =
+ g_param_spec_float ("z",
+ "Z coordinate of the reconstructed slice",
+ "Z coordinate of the reconstructed slice",
+ -G_MAXFLOAT,
+ G_MAXFLOAT,
+ 0.0f,
+ G_PARAM_READWRITE);
+
+ properties[PROP_REGION] =
+ g_param_spec_value_array ("region",
+ "Region for the parameter along z-axis as (from, to, step)",
+ "Region for the parameter along z-axis as (from, to, step)",
+ float_region_vals,
+ G_PARAM_READWRITE);
+
+ properties[PROP_PROJECTION_OFFSET] =
+ g_param_spec_value_array ("projection-offset",
+ "Offset to projection data as (x, y)",
+ "Offset to projection data as (x, y) for the case input data \
+ is cropped to the necessary range of interest",
+ region_vals,
+ G_PARAM_READWRITE);
+
+ properties[PROP_CENTER] =
+ g_param_spec_value_array ("center",
+ "Center of the volume with respect to projections (x, y)",
+ "Center of the volume with respect to projections (x, y), (rotation axes)",
+ float_region_vals,
+ G_PARAM_READWRITE);
+
+ properties[PROP_OVERALL_ANGLE] =
+ g_param_spec_float ("overall-angle",
+ "Angle covered by all projections",
+ "Angle covered by all projections (can be negative for negative steps "
+ "in case only num-projections is specified",
+ -G_MAXFLOAT,
+ G_MAXFLOAT,
+ G_PI,
+ G_PARAM_READWRITE);
+
+ properties[PROP_NUM_PROJECTIONS] =
+ g_param_spec_uint ("num-projections",
+ "Number of projections",
+ "Number of projections",
+ 0,
+ 16384,
+ 0,
+ G_PARAM_READWRITE);
+
+ properties[PROP_TOMO_ANGLE] =
+ g_param_spec_float ("tomo-angle",
+ "Tomographic rotation angle in radians",
+ "Tomographic rotation angle in radians (used for acquiring projections)",
+ -G_MAXFLOAT,
+ G_MAXFLOAT,
+ 0.0f,
+ G_PARAM_READWRITE);
+
+ properties[PROP_LAMINO_ANGLE] =
+ g_param_spec_float ("lamino-angle",
+ "Absolute laminogrpahic angle in radians",
+ "Absolute laminogrpahic angle in radians determining the sample tilt",
+ -G_MAXFLOAT,
+ G_MAXFLOAT,
+ 0.0f,
+ G_PARAM_READWRITE);
+
+ properties[PROP_ROLL_ANGLE] =
+ g_param_spec_float ("roll-angle",
+ "Sample angular misalignment to the side (roll) in radians",
+ "Sample angular misalignment to the side (roll) in radians (CW is positive)",
+ -G_MAXFLOAT,
+ G_MAXFLOAT,
+ 0.0f,
+ G_PARAM_READWRITE);
+
+ properties[PROP_PARAMETER] =
+ g_param_spec_string ("parameter",
+ "Which paramter will be varied along the z-axis",
+ "Which paramter will be varied along the z-axis, from \"z\", \"x-center\", \"lamino-angle\",\
+ \"roll-angle\"",
+ "z",
+ G_PARAM_READWRITE);
+
+ for (guint i = PROP_0 + 1; i < N_PROPERTIES; i++)
+ g_object_class_install_property (oclass, i, properties[i]);
+
+ g_type_class_add_private (oclass, sizeof(UfoLaminoBackprojectTaskPrivate));
+}
+
+static void
+ufo_lamino_backproject_task_init(UfoLaminoBackprojectTask *self)
+{
+ self->priv = UFO_LAMINO_BACKPROJECT_TASK_GET_PRIVATE(self);
+ guint i;
+ GValue int_zero = G_VALUE_INIT;
+ GValue float_zero = G_VALUE_INIT;
+
+ g_value_init (&int_zero, G_TYPE_INT);
+ g_value_init (&float_zero, G_TYPE_FLOAT);
+ g_value_set_int (&int_zero, 0);
+ g_value_set_float (&float_zero, 0.0f);
+ self->priv->x_region = g_value_array_new (3);
+ self->priv->y_region = g_value_array_new (3);
+ self->priv->region = g_value_array_new (3);
+ self->priv->z = 0.0f;
+ self->priv->projection_offset = g_value_array_new (2);
+ self->priv->center = g_value_array_new (2);
+
+ for (i = 0; i < 3; i++) {
+ g_value_array_insert (self->priv->x_region, i, &int_zero);
+ g_value_array_insert (self->priv->y_region, i, &int_zero);
+ g_value_array_insert (self->priv->region, i, &float_zero);
+ if (i < 2) {
+ g_value_array_insert (self->priv->projection_offset, i, &int_zero);
+ g_value_array_insert (self->priv->center, i, &float_zero);
+ }
+ }
+
+ self->priv->num_projections = 0;
+ self->priv->overall_angle = G_PI;
+ self->priv->tomo_angle = -G_MAXFLOAT;
+ self->priv->lamino_angle = 0.0f;
+ self->priv->roll_angle = 0.0f;
+ self->priv->parameter = PARAM_Z;
+ self->priv->count = 0;
+ self->priv->generated = FALSE;
+}
diff --git a/src/ufo-lamino-backproject-task.h b/src/ufo-lamino-backproject-task.h
new file mode 100644
index 0000000..989feaa
--- /dev/null
+++ b/src/ufo-lamino-backproject-task.h
@@ -0,0 +1,53 @@
+/*
+ * Copyright (C) 2011-2013 Karlsruhe Institute of Technology
+ *
+ * This file is part of Ufo.
+ *
+ * This library is free software: you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public
+ * License as published by the Free Software Foundation, either
+ * version 3 of the License, or (at your option) any later version.
+ *
+ * This library is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this library. If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#ifndef __UFO_LAMINO_BACKPROJECT_TASK_H
+#define __UFO_LAMINO_BACKPROJECT_TASK_H
+
+#include <ufo/ufo.h>
+
+G_BEGIN_DECLS
+
+#define UFO_TYPE_LAMINO_BACKPROJECT_TASK (ufo_lamino_backproject_task_get_type())
+#define UFO_LAMINO_BACKPROJECT_TASK(obj) (G_TYPE_CHECK_INSTANCE_CAST((obj), UFO_TYPE_LAMINO_BACKPROJECT_TASK, UfoLaminoBackprojectTask))
+#define UFO_IS_LAMINO_BACKPROJECT_TASK(obj) (G_TYPE_CHECK_INSTANCE_TYPE((obj), UFO_TYPE_LAMINO_BACKPROJECT_TASK))
+#define UFO_LAMINO_BACKPROJECT_TASK_CLASS(klass) (G_TYPE_CHECK_CLASS_CAST((klass), UFO_TYPE_LAMINO_BACKPROJECT_TASK, UfoLaminoBackprojectTaskClass))
+#define UFO_IS_LAMINO_BACKPROJECT_TASK_CLASS(klass) (G_TYPE_CHECK_CLASS_TYPE((klass), UFO_TYPE_LAMINO_BACKPROJECT_TASK))
+#define UFO_LAMINO_BACKPROJECT_TASK_GET_CLASS(obj) (G_TYPE_INSTANCE_GET_CLASS((obj), UFO_TYPE_LAMINO_BACKPROJECT_TASK, UfoLaminoBackprojectTaskClass))
+
+typedef struct _UfoLaminoBackprojectTask UfoLaminoBackprojectTask;
+typedef struct _UfoLaminoBackprojectTaskClass UfoLaminoBackprojectTaskClass;
+typedef struct _UfoLaminoBackprojectTaskPrivate UfoLaminoBackprojectTaskPrivate;
+
+struct _UfoLaminoBackprojectTask {
+ UfoTaskNode parent_instance;
+
+ UfoLaminoBackprojectTaskPrivate *priv;
+};
+
+struct _UfoLaminoBackprojectTaskClass {
+ UfoTaskNodeClass parent_class;
+};
+
+UfoNode *ufo_lamino_backproject_task_new (void);
+GType ufo_lamino_backproject_task_get_type (void);
+
+G_END_DECLS
+
+#endif \ No newline at end of file
generated by cgit v1.2.1 (git 2.18.0) at 2024-05-20 19:27:37 +0000