1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
|
/*
* Copyright (C) 2016 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 <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);
}
|
