5b6c317a
Pavel Govyadinov
implemented the c...
|
1
2
3
4
|
#ifndef STIM_CYLINDER_H
#define STIM_CYLINDER_H
#include <iostream>
#include <stim/math/circle.h>
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
5
|
#include <stim/biomodels/centerline.h>
|
e21d1051
David Mayerich
cylinder can buil...
|
6
|
#include <stim/visualization/obj.h>
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
7
|
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
8
9
10
11
|
namespace stim
{
template<typename T>
|
e21d1051
David Mayerich
cylinder can buil...
|
12
13
|
class cylinder : public centerline<T> {
protected:
|
ee17b90b
Jiaming Guo
make it compatibl...
|
14
15
|
using stim::centerline<T>::d;
|
e21d1051
David Mayerich
cylinder can buil...
|
16
|
|
e21d1051
David Mayerich
cylinder can buil...
|
17
|
std::vector< stim::vec3<T> > U; //stores the array of U vectors defining the Frenet frame
|
7115e973
David Mayerich
removed the magni...
|
18
|
std::vector< T > R; //stores a list of magnitudes for each point in the centerline (assuming mags[0] is the radius)
|
e21d1051
David Mayerich
cylinder can buil...
|
19
|
|
e21d1051
David Mayerich
cylinder can buil...
|
20
|
using stim::centerline<T>::findIdx;
|
49061f1d
Pavel Govyadinov
Modified centerli...
|
21
22
23
24
25
26
27
28
29
30
31
32
33
34
|
void init() {
U.resize(size()); //allocate space for the frenet frame vectors
// if (R.size() != 0)
R.resize(size());
stim::circle<T> c; //create a circle
stim::vec3<T> d0, d1;
for (size_t i = 0; i < size() - 1; i++) { //for each line segment in the centerline
c.rotate(d(i)); //rotate the circle to match that normal
U[i] = c.U; //save the U vector from the circle
}
U[size() - 1] = c.U; //for the last point, duplicate the final frenet frame vector
}
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
35
|
|
e21d1051
David Mayerich
cylinder can buil...
|
36
37
|
//calculates the U values for each point to initialize the frenet frame
// this function assumes that the centerline has already been set
|
e21d1051
David Mayerich
cylinder can buil...
|
38
|
|
e21d1051
David Mayerich
cylinder can buil...
|
39
|
public:
|
e21d1051
David Mayerich
cylinder can buil...
|
40
|
|
eb554b48
David Mayerich
bug fixes for Net...
|
41
42
|
using stim::centerline<T>::size;
using stim::centerline<T>::at;
|
ee17b90b
Jiaming Guo
make it compatibl...
|
43
44
|
using stim::centerline<T>::L;
using stim::centerline<T>::length;
|
eb554b48
David Mayerich
bug fixes for Net...
|
45
|
|
ee17b90b
Jiaming Guo
make it compatibl...
|
46
|
cylinder() : centerline<T>(){}
|
e21d1051
David Mayerich
cylinder can buil...
|
47
|
|
ee17b90b
Jiaming Guo
make it compatibl...
|
48
|
cylinder(centerline<T>c) : centerline<T>(c) {
|
e21d1051
David Mayerich
cylinder can buil...
|
49
|
init();
|
e21d1051
David Mayerich
cylinder can buil...
|
50
|
}
|
49061f1d
Pavel Govyadinov
Modified centerli...
|
51
52
53
54
55
56
|
cylinder(std::vector<stim::vec3<T> > p, std::vector<T> s)
: centerline<T>(p)
{
R = s;
init();
|
6aa276f6
Pavel Govyadinov
Parentheses fixed
|
57
|
}
|
49061f1d
Pavel Govyadinov
Modified centerli...
|
58
59
60
|
cylinder(stim::centerline<T> p, std::vector<T> s)
{
|
329587ef
David Mayerich
fixed compiler er...
|
61
62
|
//was d = s;
p = s;
|
49061f1d
Pavel Govyadinov
Modified centerli...
|
63
|
init();
|
6aa276f6
Pavel Govyadinov
Parentheses fixed
|
64
|
}
|
3b14b0f5
Pavel Govyadinov
added operators f...
|
65
|
|
ee17b90b
Jiaming Guo
make it compatibl...
|
66
67
68
69
70
|
//cylinder(centerline<T>c, T r) : centerline(c) {
// init();
// //add_mag(r);
//}
|
e21d1051
David Mayerich
cylinder can buil...
|
71
|
//initialize a cylinder with a list of points and magnitude values
|
ee17b90b
Jiaming Guo
make it compatibl...
|
72
73
74
75
|
//cylinder(centerline<T>c, std::vector<T> r) : centerline(c){
// init();
// //add_mag(r);
//}
|
e21d1051
David Mayerich
cylinder can buil...
|
76
|
|
58674036
Jiaming Guo
fix minor error f...
|
77
78
79
80
81
82
|
//copy the original radius
void copy_r(std::vector<T> radius) {
for (unsigned i = 0; i < radius.size(); i++)
R[i] = radius[i];
}
|
e21d1051
David Mayerich
cylinder can buil...
|
83
84
85
86
|
///Returns magnitude i at the given length into the fiber (based on the pvalue).
///Interpolates the position along the line.
///@param l: the location of the in the cylinder.
///@param idx: integer location of the point closest to l but prior to it.
|
7115e973
David Mayerich
removed the magni...
|
87
|
T r(T l, int idx) {
|
e21d1051
David Mayerich
cylinder can buil...
|
88
|
T a = (l - L[idx]) / (L[idx + 1] - L[idx]);
|
7115e973
David Mayerich
removed the magni...
|
89
|
T v2 = (R[idx] + (R[idx + 1] - R[idx])*a);
|
e21d1051
David Mayerich
cylinder can buil...
|
90
91
92
93
94
95
96
|
return v2;
}
///Returns the ith magnitude at the given p-value (p value ranges from 0 to 1).
///interpolates the radius along the line.
///@param pvalue: the location of the in the cylinder, from 0 (beginning to 1).
|
720240c8
Jiaming Guo
everything now is...
|
97
|
T rl(T pvalue) {
|
7115e973
David Mayerich
removed the magni...
|
98
99
|
if (pvalue <= 0.0) return R[0];
if (pvalue >= 1.0) return R[size() - 1];
|
e21d1051
David Mayerich
cylinder can buil...
|
100
101
102
|
T l = pvalue*L[L.size() - 1];
int idx = findIdx(l);
|
0d0ef1a9
David Mayerich
Adapted classes t...
|
103
|
return r(l, idx);
|
e21d1051
David Mayerich
cylinder can buil...
|
104
105
|
}
|
eb554b48
David Mayerich
bug fixes for Net...
|
106
107
|
/// Returns the magnitude at the given index
/// @param i is the index of the desired point
|
7115e973
David Mayerich
removed the magni...
|
108
|
/// @param r is the index of the magnitude value
|
720240c8
Jiaming Guo
everything now is...
|
109
110
111
|
T r(unsigned i) {
return R[i];
}
|
eb554b48
David Mayerich
bug fixes for Net...
|
112
|
|
0e0feff0
Jiaming Guo
fixed bugs in mer...
|
113
|
|
eb554b48
David Mayerich
bug fixes for Net...
|
114
|
///adds a magnitude to each point in the cylinder
|
7115e973
David Mayerich
removed the magni...
|
115
|
/*void add_mag(V val = 0) {
|
eb554b48
David Mayerich
bug fixes for Net...
|
116
117
|
if (M.size() == 0) M.resize(size()); //if the magnitude vector isn't initialized, resize it to match the centerline
for (size_t i = 0; i < size(); i++) //for each point
|
7115e973
David Mayerich
removed the magni...
|
118
119
|
R[i].push_back(val); //add this value to the magnitude vector at each point
}*/
|
eb554b48
David Mayerich
bug fixes for Net...
|
120
121
|
//adds a magnitude based on a list of magnitudes for each point
|
7115e973
David Mayerich
removed the magni...
|
122
|
/*void add_mag(std::vector<T> val) {
|
eb554b48
David Mayerich
bug fixes for Net...
|
123
124
|
if (M.size() == 0) M.resize(size()); //if the magnitude vector isn't initialized, resize it to match the centerline
for (size_t i = 0; i < size(); i++) //for each point
|
7115e973
David Mayerich
removed the magni...
|
125
126
|
R[i].push_back(val[i]); //add this value to the magnitude vector at each point
}*/
|
eb554b48
David Mayerich
bug fixes for Net...
|
127
128
|
//sets the value of magnitude m at point i
|
7115e973
David Mayerich
removed the magni...
|
129
130
|
void set_r(size_t i, T r) {
R[i] = r;
|
eb554b48
David Mayerich
bug fixes for Net...
|
131
132
|
}
|
7115e973
David Mayerich
removed the magni...
|
133
|
/*size_t nmags() {
|
eb554b48
David Mayerich
bug fixes for Net...
|
134
|
if (M.size() == 0) return 0;
|
7115e973
David Mayerich
removed the magni...
|
135
136
|
else return R[0].size();
}*/
|
eb554b48
David Mayerich
bug fixes for Net...
|
137
|
|
e21d1051
David Mayerich
cylinder can buil...
|
138
|
///Returns a circle representing the cylinder cross section at point i
|
7115e973
David Mayerich
removed the magni...
|
139
140
|
stim::circle<T> circ(size_t i) {
return stim::circle<T>(at(i), R[i], d(i), U[i]);
|
e21d1051
David Mayerich
cylinder can buil...
|
141
142
143
|
}
///Returns an OBJ object representing the cylinder with a radial tesselation value of N using magnitude m
|
7115e973
David Mayerich
removed the magni...
|
144
|
stim::obj<T> OBJ(size_t N) {
|
e21d1051
David Mayerich
cylinder can buil...
|
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
|
stim::obj<T> out; //create an OBJ object
stim::circle<T> c0, c1;
std::vector< stim::vec3<T> > p0, p1; //allocate space for the point sets representing the circles bounding each cylinder segment
T u0, u1, v0, v1; //allocate variables to store running texture coordinates
for (size_t i = 1; i < size(); i++) { //for each line segment in the cylinder
c0 = circ(i - 1); //get the two circles bounding the line segment
c1 = circ(i);
p0 = c0.points(N); //get t points for each of the end caps
p1 = c1.points(N);
u0 = L[i - 1] / length(); //calculate the texture coordinate (u, v) where u runs along the cylinder length
u1 = L[i] / length();
for (size_t n = 1; n < N; n++) { //for each point in the circle
v0 = (double)(n-1) / (double)(N - 1); //v texture coordinate runs around the cylinder
v1 = (double)(n) / (double)(N - 1);
out.Begin(OBJ_FACE); //start a face (quad)
out.TexCoord(u0, v0);
out.Vertex(p0[n - 1][0], p0[n - 1][1], p0[n - 1][2]); //output the points composing a strip of quads wrapping the cylinder segment
out.TexCoord(u1, v0);
out.Vertex(p1[n - 1][0], p1[n - 1][1], p1[n - 1][2]);
out.TexCoord(u0, v1);
out.Vertex(p1[n][0], p1[n][1], p1[n][2]);
out.TexCoord(u1, v1);
out.Vertex(p0[n][0], p0[n][1], p0[n][2]);
out.End();
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
173
|
}
|
e21d1051
David Mayerich
cylinder can buil...
|
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
|
v0 = (double)(N - 2) / (double)(N - 1); //v texture coordinate runs around the cylinder
v1 = 1.0;
out.Begin(OBJ_FACE);
out.TexCoord(u0, v0);
out.Vertex(p0[N - 1][0], p0[N - 1][1], p0[N - 1][2]); //output the points composing a strip of quads wrapping the cylinder segment
out.TexCoord(u1, v0);
out.Vertex(p1[N - 1][0], p1[N - 1][1], p1[N - 1][2]);
out.TexCoord(u0, v1);
out.Vertex(p1[0][0], p1[0][1], p1[0][2]);
out.TexCoord(u1, v1);
out.Vertex(p0[0][0], p0[0][1], p0[0][2]);
out.End();
}
return out;
}
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
190
|
|
eb554b48
David Mayerich
bug fixes for Net...
|
191
192
193
194
195
|
std::string str() {
std::stringstream ss;
size_t N = std::vector< stim::vec3<T> >::size();
ss << "---------[" << N << "]---------" << std::endl;
for (size_t i = 0; i < N; i++)
|
7115e973
David Mayerich
removed the magni...
|
196
|
ss << std::vector< stim::vec3<T> >::at(i) << " r = " << R[i] << " u = " << U[i] << std::endl;
|
eb554b48
David Mayerich
bug fixes for Net...
|
197
198
199
200
201
202
203
|
ss << "--------------------" << std::endl;
return ss.str();
}
/// Integrates a magnitude value along the cylinder.
/// @param m is the magnitude value to be integrated (this is usually the radius)
|
7115e973
David Mayerich
removed the magni...
|
204
|
T integrate() {
|
58674036
Jiaming Guo
fix minor error f...
|
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
|
T sum = 0; //initialize the integral to zero
if (L.size() == 1)
return sum;
else {
T m0, m1; //allocate space for both magnitudes in a single segment
m0 = R[0]; //initialize the first point and magnitude to the first point in the cylinder
T len = L[1]; //allocate space for the segment length
for (unsigned p = 1; p < size(); p++) { //for every consecutive point in the cylinder
m1 = R[p];
if (p > 1) len = (L[p] - L[p - 1]); //calculate the segment length using the L array
sum += (m0 + m1) / (T)2.0 * len; //add the average magnitude, weighted by the segment length
m0 = m1; //move to the next segment by shifting points
}
return sum; //return the integral
|
eb554b48
David Mayerich
bug fixes for Net...
|
221
|
}
|
eb554b48
David Mayerich
bug fixes for Net...
|
222
223
224
225
226
227
228
229
|
}
/// Resamples the cylinder to provide a maximum distance of "spacing" between centerline points. All current
/// centerline points are guaranteed to exist in the new cylinder
/// @param spacing is the maximum spacing allowed between sample points
cylinder<T> resample(T spacing) {
cylinder<T> c = stim::centerline<T>::resample(spacing); //resample the centerline and use it to create a new cylinder
|
7115e973
David Mayerich
removed the magni...
|
230
231
232
233
234
|
//size_t nm = nmags(); //get the number of magnitude values in the current cylinder
//if (nm > 0) { //if there are magnitude values
// std::vector<T> magvec(nm, 0); //create a magnitude vector for a single point
// c.M.resize(c.size(), magvec); //allocate space for a magnitude vector at each point of the new cylinder
//}
|
eb554b48
David Mayerich
bug fixes for Net...
|
235
236
237
238
239
|
T l, t;
for (size_t i = 0; i < c.size(); i++) { //for each point in the new cylinder
l = c.L[i]; //get the length along the new cylinder
t = l / length(); //calculate the parameter value along the new cylinder
|
7115e973
David Mayerich
removed the magni...
|
240
241
242
|
//for (size_t mag = 0; mag < nm; mag++) { //for each magnitude value
c.R[i] = r(t); //retrieve the interpolated magnitude from the current cylinder and store it in the new one
//}
|
eb554b48
David Mayerich
bug fixes for Net...
|
243
|
}
|
0d0ef1a9
David Mayerich
Adapted classes t...
|
244
245
|
return c;
}
|
eb554b48
David Mayerich
bug fixes for Net...
|
246
|
|
0d0ef1a9
David Mayerich
Adapted classes t...
|
247
248
249
250
251
252
253
254
255
256
|
std::vector< stim::cylinder<T> > split(unsigned int idx) {
unsigned N = size();
std::vector< stim::centerline<T> > LL;
LL.resize(2);
LL = (*this).centerline<T>::split(idx);
std::vector< stim::cylinder<T> > C(LL.size());
unsigned i = 0;
C[0] = LL[0];
//C[0].R.resize(idx);
|
da89bce9
Jiaming Guo
fixed splitting bugs
|
257
|
for (; i < idx + 1; i++) {
|
0d0ef1a9
David Mayerich
Adapted classes t...
|
258
259
260
261
262
263
264
|
//for(unsigned d = 0; d < 3; d++)
//C[0][i][d] = LL[0].c[i][d];
C[0].R[i] = R[i];
//C[0].R[i].resize(1);
}
if (C.size() == 2) {
C[1] = LL[1];
|
da89bce9
Jiaming Guo
fixed splitting bugs
|
265
|
i--;
|
0d0ef1a9
David Mayerich
Adapted classes t...
|
266
267
268
269
270
271
272
273
|
//C[1].M.resize(N - idx);
for (; i < N; i++) {
//for(unsigned d = 0; d < 3; d++)
//C[1][i - idx][d] = LL[1].c[i - idx][d];
C[1].R[i - idx] = R[i];
//C[1].M[i - idx].resize(1);
}
}
|
eb554b48
David Mayerich
bug fixes for Net...
|
274
|
|
0d0ef1a9
David Mayerich
Adapted classes t...
|
275
|
return C;
|
eb554b48
David Mayerich
bug fixes for Net...
|
276
277
|
}
|
e21d1051
David Mayerich
cylinder can buil...
|
278
279
280
|
/*
///inits the cylinder from a list of points (std::vector of stim::vec3 --inP) and magnitudes (inM)
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
281
|
void
|
e21d1051
David Mayerich
cylinder can buil...
|
282
283
284
|
init(centerline inP, std::vector< std::vector<T> > inM) {
M = inM; //the magnitude vector can be copied directly
(*this) = inP; //the centerline can be copied to this class directly
|
308a743c
David Mayerich
fixed class compa...
|
285
286
|
stim::vec3<float> v1;
stim::vec3<float> v2;
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
287
|
e.resize(inP.size());
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
288
289
290
291
|
norms.resize(inP.size());
Us.resize(inP.size());
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
292
293
|
if(inP.size() < 2)
return;
|
8495a970
Pavel Govyadinov
added comments to...
|
294
295
|
//calculate each L.
|
1b2858a2
David Mayerich
fixed some commen...
|
296
297
298
|
L.resize(inP.size()); //the number of precomputed lengths will equal the number of points
T temp = (T)0; //length up to that point
L[0] = temp;
|
308a743c
David Mayerich
fixed class compa...
|
299
|
for(size_t i = 1; i < L.size(); i++)
|
10654a1f
Pavel Govyadinov
added the necessa...
|
300
|
{
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
301
|
temp += (inP[i-1] - inP[i]).len();
|
bf23ee36
Pavel Govyadinov
more minor bug fi...
|
302
|
L[i] = temp;
|
10654a1f
Pavel Govyadinov
added the necessa...
|
303
|
}
|
9c97e126
David Mayerich
added an axis-ali...
|
304
|
|
308a743c
David Mayerich
fixed class compa...
|
305
|
stim::vec3<T> dr = (inP[1] - inP[0]).norm();
|
7115e973
David Mayerich
removed the magni...
|
306
|
s = stim::circle<T>(inP[0], inR[0][0], dr, stim::vec3<T>(1,0,0));
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
307
|
norms[0] = s.N;
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
308
|
e[0] = s;
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
309
|
Us[0] = e[0].U;
|
308a743c
David Mayerich
fixed class compa...
|
310
|
for(size_t i = 1; i < inP.size()-1; i++)
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
311
312
313
314
315
316
|
{
s.center(inP[i]);
v1 = (inP[i] - inP[i-1]).norm();
v2 = (inP[i+1] - inP[i]).norm();
dr = (v1+v2).norm();
s.normal(dr);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
317
318
319
|
norms[i] = s.N;
|
7115e973
David Mayerich
removed the magni...
|
320
|
s.scale(inR[i][0]/inR[i-1][0]);
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
321
|
e[i] = s;
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
322
|
Us[i] = e[i].U;
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
323
324
325
326
327
328
|
}
int j = inP.size()-1;
s.center(inP[j]);
dr = (inP[j] - inP[j-1]).norm();
s.normal(dr);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
329
330
331
|
norms[j] = s.N;
|
7115e973
David Mayerich
removed the magni...
|
332
|
s.scale(inR[j][0]/inR[j-1][0]);
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
333
|
e[j] = s;
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
334
|
Us[j] = e[j].U;
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
335
336
|
}
|
8495a970
Pavel Govyadinov
added comments to...
|
337
|
///returns the direction vector at point idx.
|
308a743c
David Mayerich
fixed class compa...
|
338
|
stim::vec3<T>
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
339
340
341
342
|
d(int idx)
{
if(idx == 0)
{
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
343
344
345
346
347
348
349
|
stim::vec3<T> temp(
c[idx+1][0]-c[idx][0],
c[idx+1][1]-c[idx][1],
c[idx+1][2]-c[idx][2]
);
// return (e[idx+1].P - e[idx].P).norm();
return (temp.norm());
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
350
|
}
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
351
|
else if(idx == N-1)
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
352
|
{
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
353
354
355
356
357
358
359
|
stim::vec3<T> temp(
c[idx][0]-c[idx+1][0],
c[idx][1]-c[idx+1][1],
c[idx][2]-c[idx+1][2]
);
// return (e[idx].P - e[idx-1].P).norm();
return (temp.norm());
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
360
361
362
|
}
else
{
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
363
|
// return (e[idx+1].P - e[idx].P).norm();
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
|
// stim::vec3<float> v1 = (e[idx].P-e[idx-1].P).norm();
stim::vec3<T> v1(
c[idx][0]-c[idx-1][0],
c[idx][1]-c[idx-1][1],
c[idx][2]-c[idx-1][2]
);
// stim::vec3<float> v2 = (e[idx+1].P-e[idx].P).norm();
stim::vec3<T> v2(
c[idx+1][0]-c[idx][0],
c[idx+1][1]-c[idx][1],
c[idx+1][2]-c[idx][2]
);
return (v1.norm()+v2.norm()).norm();
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
379
380
381
|
}
// return e[idx].N;
|
10654a1f
Pavel Govyadinov
added the necessa...
|
382
383
|
}
|
308a743c
David Mayerich
fixed class compa...
|
384
|
stim::vec3<T>
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
385
|
d(T l, int idx)
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
386
|
{
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
387
|
if(idx == 0 || idx == N-1)
|
10654a1f
Pavel Govyadinov
added the necessa...
|
388
|
{
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
389
390
|
return norms[idx];
// return e[idx].N;
|
10654a1f
Pavel Govyadinov
added the necessa...
|
391
|
}
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
392
393
|
else
{
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
394
|
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
395
|
T rat = (l-L[idx])/(L[idx+1]-L[idx]);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
396
397
|
return( norms[idx] + (norms[idx+1] - norms[idx])*rat);
// return( e[idx].N + (e[idx+1].N - e[idx].N)*rat);
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
398
|
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
399
400
|
}
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
401
|
|
8495a970
Pavel Govyadinov
added comments to...
|
402
403
|
///finds the index of the point closest to the length l on the lower bound.
///binary search.
|
10654a1f
Pavel Govyadinov
added the necessa...
|
404
|
int
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
405
406
407
408
409
|
findIdx(T l)
{
unsigned int i = L.size()/2;
unsigned int max = L.size()-1;
unsigned int min = 0;
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
410
|
while(i > 0 && i < L.size()-1)
|
10654a1f
Pavel Govyadinov
added the necessa...
|
411
|
{
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
412
413
|
// std::cerr << "Trying " << i << std::endl;
// std::cerr << "l is " << l << ", L[" << i << "]" << L[i] << std::endl;
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
414
|
if(l < L[i])
|
10654a1f
Pavel Govyadinov
added the necessa...
|
415
|
{
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
416
|
max = i;
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
417
|
i = min+(max-min)/2;
|
10654a1f
Pavel Govyadinov
added the necessa...
|
418
|
}
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
419
|
else if(L[i] <= l && L[i+1] >= l)
|
10654a1f
Pavel Govyadinov
added the necessa...
|
420
421
422
423
424
|
{
break;
}
else
{
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
425
|
min = i;
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
426
|
i = min+(max-min)/2;
|
10654a1f
Pavel Govyadinov
added the necessa...
|
427
428
429
430
|
}
}
return i;
}
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
431
|
|
9c97e126
David Mayerich
added an axis-ali...
|
432
433
|
public:
///default constructor
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
434
|
cylinder()
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
435
|
// : centerline<T>()
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
436
437
438
|
{
}
|
9c97e126
David Mayerich
added an axis-ali...
|
439
|
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
440
|
///constructor to create a cylinder from a set of points, radii, and the number of sides for the cylinder.
|
7d3162a2
Pavel Govyadinov
fixed majority of...
|
441
|
///@param inP: Vector of stim vec3 composing the points of the centerline.
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
442
|
///@param inM: Vector of stim vecs composing the radii of the centerline.
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
443
444
445
|
cylinder(std::vector<stim::vec3<T> > inP, std::vector<stim::vec<T> > inM)
: centerline<T>(inP)
{
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
446
447
448
|
init(inP, inM);
}
|
09049866
Pavel Govyadinov
fixed a drawing b...
|
449
450
451
452
453
454
455
|
///constructor to create a cylinder from a set of points, radii, and the number of sides for the cylinder.
///@param inP: Vector of stim vec3 composing the points of the centerline.
///@param inM: Vector of stim vecs composing the radii of the centerline.
cylinder(std::vector<stim::vec3<T> > inP, std::vector< T > inM)
: centerline<T>(inP)
{
std::vector<stim::vec<T> > temp;
|
a36c3850
Pavel Govyadinov
fixed the bug tha...
|
456
|
stim::vec<T> zero(0.0,0.0);
|
09049866
Pavel Govyadinov
fixed a drawing b...
|
457
|
temp.resize(inM.size(), zero);
|
47fe6300
Pavel Govyadinov
fixed an accident...
|
458
|
for(int i = 0; i < inM.size(); i++)
|
7115e973
David Mayerich
removed the magni...
|
459
|
temp[i][0] = inR[i];
|
09049866
Pavel Govyadinov
fixed a drawing b...
|
460
461
462
|
init(inP, temp);
}
|
7d3162a2
Pavel Govyadinov
fixed majority of...
|
463
|
|
58ee2b21
David Mayerich
incorporated stim...
|
464
|
///Constructor defines a cylinder with centerline inP and magnitudes of zero
|
7d3162a2
Pavel Govyadinov
fixed majority of...
|
465
|
///@param inP: Vector of stim vec3 composing the points of the centerline
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
466
467
468
|
cylinder(std::vector< stim::vec3<T> > inP)
: centerline<T>(inP)
{
|
96e8ab1f
David Mayerich
fixed stim::cylin...
|
469
|
std::vector< stim::vec<T> > inM; //create an array of arbitrary magnitudes
|
9c97e126
David Mayerich
added an axis-ali...
|
470
471
472
473
474
|
stim::vec<T> zero;
zero.push_back(0);
inM.resize(inP.size(), zero); //initialize the magnitude values to zero
|
58ee2b21
David Mayerich
incorporated stim...
|
475
476
477
|
init(inP, inM);
}
|
e21d1051
David Mayerich
cylinder can buil...
|
478
479
480
481
482
|
//assignment operator creates a cylinder from a centerline (default radius is zero)
cylinder& operator=(stim::centerline<T> c) {
init(c);
}
|
58ee2b21
David Mayerich
incorporated stim...
|
483
484
485
|
///Returns the number of points on the cylinder centerline
unsigned int size(){
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
486
|
return N;
|
58ee2b21
David Mayerich
incorporated stim...
|
487
488
|
}
|
e21d1051
David Mayerich
cylinder can buil...
|
489
|
|
10654a1f
Pavel Govyadinov
added the necessa...
|
490
|
///Returns a position vector at the given p-value (p value ranges from 0 to 1).
|
8495a970
Pavel Govyadinov
added comments to...
|
491
492
|
///interpolates the position along the line.
///@param pvalue: the location of the in the cylinder, from 0 (beginning to 1).
|
308a743c
David Mayerich
fixed class compa...
|
493
|
stim::vec3<T>
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
494
495
|
p(T pvalue)
{
|
bf23ee36
Pavel Govyadinov
more minor bug fi...
|
496
|
if(pvalue < 0.0 || pvalue > 1.0)
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
497
|
{
|
308a743c
David Mayerich
fixed class compa...
|
498
|
return stim::vec3<float>(-1,-1,-1);
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
499
|
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
500
501
|
T l = pvalue*L[L.size()-1];
int idx = findIdx(l);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
502
|
return (p(l,idx));
|
e21d1051
David Mayerich
cylinder can buil...
|
503
|
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
504
|
|
8495a970
Pavel Govyadinov
added comments to...
|
505
506
507
508
|
///Returns a position vector at the given length into the fiber (based on the pvalue).
///Interpolates the radius along the line.
///@param l: the location of the in the cylinder.
///@param idx: integer location of the point closest to l but prior to it.
|
308a743c
David Mayerich
fixed class compa...
|
509
|
stim::vec3<T>
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
510
511
512
|
p(T l, int idx)
{
T rat = (l-L[idx])/(L[idx+1]-L[idx]);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
513
514
515
516
517
518
519
520
521
522
523
524
525
|
stim::vec3<T> v1(
c[idx][0],
c[idx][1],
c[idx][2]
);
stim::vec3<T> v2(
c[idx+1][0],
c[idx+1][1],
c[idx+1][2]
);
// return( e[idx].P + (e[idx+1].P-e[idx].P)*rat);
return( v1 + (v2-v1)*rat);
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
526
527
|
// return(
// return (pos[idx] + (pos[idx+1]-pos[idx])*((l-L[idx])/(L[idx+1]- L[idx])));
|
10654a1f
Pavel Govyadinov
added the necessa...
|
528
529
530
|
}
///Returns a radius at the given p-value (p value ranges from 0 to 1).
|
8495a970
Pavel Govyadinov
added comments to...
|
531
532
|
///interpolates the radius along the line.
///@param pvalue: the location of the in the cylinder, from 0 (beginning to 1).
|
10654a1f
Pavel Govyadinov
added the necessa...
|
533
|
T
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
534
535
|
r(T pvalue)
{
|
98eecaa9
David Mayerich
VS and win32 updates
|
536
537
538
539
|
if(pvalue < 0.0 || pvalue > 1.0){
std::cerr<<"Error, value "<<pvalue<<" is outside of [0 1]."<<std::endl;
exit(1);
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
540
541
|
T l = pvalue*L[L.size()-1];
int idx = findIdx(l);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
542
|
return (r(l,idx));
|
e21d1051
David Mayerich
cylinder can buil...
|
543
|
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
544
|
|
8495a970
Pavel Govyadinov
added comments to...
|
545
546
547
548
|
///Returns a radius at the given length into the fiber (based on the pvalue).
///Interpolates the position along the line.
///@param l: the location of the in the cylinder.
///@param idx: integer location of the point closest to l but prior to it.
|
10654a1f
Pavel Govyadinov
added the necessa...
|
549
|
T
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
550
551
552
|
r(T l, int idx)
{
T rat = (l-L[idx])/(L[idx+1]-L[idx]);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
553
554
555
556
557
558
559
560
|
T v1 = (e[idx].U.len() + (e[idx+1].U.len() - e[idx].U.len())*rat);
T v3 = (Us[idx].len() + (Us[idx+1].len() - Us[idx].len())*rat);
T v2 = (mags[idx][0] + (mags[idx+1][0]-mags[idx][0])*rat);
// std::cout << (float)v1 = (float) v2 << std::endl;
if(abs(v3 - v1) >= 10e-6)
{
std::cout << "-------------------------" << std::endl;
std::cout << e[idx].str() << std::endl << std::endl;
|
b50c840e
David Mayerich
eliminated ANN fr...
|
561
|
std::cout << Us[idx].str() << std::endl;
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
562
563
564
565
566
|
std::cout << (float)v1 - (float) v2 << std::endl;
std::cout << "failed" << std::endl;
}
// std::cout << e[idx].U.len() << " " << mags[idx][0] << std::endl;
// std::cout << v2 << std::endl;
|
24dee661
Pavel Govyadinov
finished the debu...
|
567
|
return(v2);
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
568
569
|
// return (mags[idx][0] + (mags[idx+1][0]-mags[idx][0])*rat);
// (
|
66fe63f3
Pavel Govyadinov
fixed minor typo
|
570
|
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
571
|
|
9c97e126
David Mayerich
added an axis-ali...
|
572
573
574
575
576
577
578
579
580
581
|
/// Returns the magnitude at the given index
/// @param i is the index of the desired point
/// @param m is the index of the magnitude value
T ri(unsigned i, unsigned m = 0){
return mags[i][m];
}
/// Adds a new magnitude value to all points
/// @param m is the starting value for the new magnitude
void add_mag(T m = 0){
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
582
|
for(unsigned int p = 0; p < N; p++)
|
9c97e126
David Mayerich
added an axis-ali...
|
583
584
585
586
587
588
589
590
591
592
593
|
mags[p].push_back(m);
}
/// Sets a magnitude value
/// @param val is the new value for the magnitude
/// @param p is the point index for the magnitude to be set
/// @param m is the index for the magnitude
void set_mag(T val, unsigned p, unsigned m = 0){
mags[p][m] = val;
}
|
b3a38641
David Mayerich
added the ability...
|
594
595
596
597
598
|
/// Returns the number of magnitude values at each point
unsigned nmags(){
return mags[0].size();
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
599
|
///returns the position of the point with a given pvalue and theta on the surface
|
8495a970
Pavel Govyadinov
added comments to...
|
600
601
602
|
///in x, y, z coordinates. Theta is in degrees from 0 to 360.
///@param pvalue: the location of the in the cylinder, from 0 (beginning to 1).
///@param theta: the angle to the point of a circle.
|
308a743c
David Mayerich
fixed class compa...
|
603
|
stim::vec3<T>
|
10654a1f
Pavel Govyadinov
added the necessa...
|
604
605
|
surf(T pvalue, T theta)
{
|
bf23ee36
Pavel Govyadinov
more minor bug fi...
|
606
|
if(pvalue < 0.0 || pvalue > 1.0)
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
607
|
{
|
308a743c
David Mayerich
fixed class compa...
|
608
|
return stim::vec3<float>(-1,-1,-1);
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
609
|
} else {
|
10654a1f
Pavel Govyadinov
added the necessa...
|
610
611
|
T l = pvalue*L[L.size()-1];
int idx = findIdx(l);
|
308a743c
David Mayerich
fixed class compa...
|
612
|
stim::vec3<T> ps = p(l, idx);
|
10654a1f
Pavel Govyadinov
added the necessa...
|
613
|
T m = r(l, idx);
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
614
615
616
617
|
s = e[idx];
s.center(ps);
s.normal(d(l, idx));
s.scale(m/e[idx].U.len());
|
10654a1f
Pavel Govyadinov
added the necessa...
|
618
|
return(s.p(theta));
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
619
|
}
|
10654a1f
Pavel Govyadinov
added the necessa...
|
620
621
|
}
|
8495a970
Pavel Govyadinov
added comments to...
|
622
|
///returns a vector of points necessary to create a circle at every position in the fiber.
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
623
|
///@param sides: the number of sides of each circle.
|
308a743c
David Mayerich
fixed class compa...
|
624
|
std::vector<std::vector<vec3<T> > >
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
625
626
|
getPoints(int sides)
{
|
308a743c
David Mayerich
fixed class compa...
|
627
|
std::vector<std::vector <vec3<T> > > points;
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
628
629
|
points.resize(N);
for(int i = 0; i < N; i++)
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
630
|
{
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
631
|
points[i] = e[i].getPoints(sides);
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
632
|
}
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
633
|
return points;
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
634
|
}
|
58ee2b21
David Mayerich
incorporated stim...
|
635
|
|
8c4f5d84
Pavel Govyadinov
fixed the issue w...
|
636
637
638
639
640
641
|
///returns the total length of the line at index j.
T
getl(int j)
{
return (L[j]);
}
|
58ee2b21
David Mayerich
incorporated stim...
|
642
643
|
/// Allows a point on the centerline to be accessed using bracket notation
|
308a743c
David Mayerich
fixed class compa...
|
644
|
vec3<T> operator[](unsigned int i){
|
3e5d3ad3
Pavel Govyadinov
merged the change...
|
645
|
return e[i].P;
|
58ee2b21
David Mayerich
incorporated stim...
|
646
647
648
649
650
651
652
|
}
/// Returns the total length of the cylinder centerline
T length(){
return L.back();
}
|
9c97e126
David Mayerich
added an axis-ali...
|
653
654
655
656
657
658
659
|
/// Integrates a magnitude value along the cylinder.
/// @param m is the magnitude value to be integrated (this is usually the radius)
T integrate(unsigned m = 0){
T M = 0; //initialize the integral to zero
T m0, m1; //allocate space for both magnitudes in a single segment
|
308a743c
David Mayerich
fixed class compa...
|
660
|
//vec3<T> p0, p1; //allocate space for both points in a single segment
|
9c97e126
David Mayerich
added an axis-ali...
|
661
662
663
664
665
666
667
|
m0 = mags[0][m]; //initialize the first point and magnitude to the first point in the cylinder
//p0 = pos[0];
T len = L[0]; //allocate space for the segment length
//for every consecutive point in the cylinder
|
0311ab81
Pavel Govyadinov
fixed some issues...
|
668
|
for(unsigned p = 1; p < N; p++){
|
9c97e126
David Mayerich
added an axis-ali...
|
669
670
671
672
673
674
675
|
//p1 = pos[p]; //get the position and magnitude for the next point
m1 = mags[p][m];
if(p > 1) len = (L[p-1] - L[p-2]); //calculate the segment length using the L array
//add the average magnitude, weighted by the segment length
|
308a743c
David Mayerich
fixed class compa...
|
676
|
M += (m0 + m1)/(T)2.0 * len;
|
9c97e126
David Mayerich
added an axis-ali...
|
677
678
679
680
681
682
683
684
685
686
687
688
689
690
|
m0 = m1; //move to the next segment by shifting points
}
return M; //return the integral
}
/// Averages a magnitude value across the cylinder
/// @param m is the magnitude value to be averaged (this is usually the radius)
T average(unsigned m = 0){
//return the average magnitude
return integrate(m) / L.back();
}
|
58ee2b21
David Mayerich
incorporated stim...
|
691
692
693
694
695
|
/// Resamples the cylinder to provide a maximum distance of "spacing" between centerline points. All current
/// centerline points are guaranteed to exist in the new cylinder
/// @param spacing is the maximum spacing allowed between sample points
cylinder<T> resample(T spacing){
|
308a743c
David Mayerich
fixed class compa...
|
696
|
std::vector< vec3<T> > result;
|
9c97e126
David Mayerich
added an axis-ali...
|
697
|
|
1b2858a2
David Mayerich
fixed some commen...
|
698
|
vec3<T> p0 = e[0].P; //initialize p0 to the first point on the centerline
|
308a743c
David Mayerich
fixed class compa...
|
699
|
vec3<T> p1;
|
1b2858a2
David Mayerich
fixed some commen...
|
700
|
unsigned N = size(); //number of points in the current centerline
|
58ee2b21
David Mayerich
incorporated stim...
|
701
702
703
|
//for each line segment on the centerline
for(unsigned int i = 1; i < N; i++){
|
1b2858a2
David Mayerich
fixed some commen...
|
704
|
p1 = e[i].P; //get the second point in the line segment
|
58ee2b21
David Mayerich
incorporated stim...
|
705
|
|
1b2858a2
David Mayerich
fixed some commen...
|
706
707
|
vec3<T> v = p1 - p0; //calculate the vector between these two points
T d = v.len(); //calculate the distance between these two points (length of the line segment)
|
58ee2b21
David Mayerich
incorporated stim...
|
708
|
|
308a743c
David Mayerich
fixed class compa...
|
709
|
size_t nsteps = (size_t)std::ceil(d / spacing); //calculate the number of steps to take along the segment to meet the spacing criteria
|
1b2858a2
David Mayerich
fixed some commen...
|
710
|
T stepsize = (T)1.0 / nsteps; //calculate the parametric step size between new centerline points
|
58ee2b21
David Mayerich
incorporated stim...
|
711
712
713
|
//for each step along the line segment
for(unsigned s = 0; s < nsteps; s++){
|
1b2858a2
David Mayerich
fixed some commen...
|
714
715
|
T alpha = stepsize * s; //calculate the fraction of the distance along the line segment covered
result.push_back(p0 + alpha * v); //push the point at alpha position along the line segment
|
58ee2b21
David Mayerich
incorporated stim...
|
716
717
|
}
|
1b2858a2
David Mayerich
fixed some commen...
|
718
|
p0 = p1; //shift the points to move to the next line segment
|
58ee2b21
David Mayerich
incorporated stim...
|
719
720
|
}
|
1b2858a2
David Mayerich
fixed some commen...
|
721
|
result.push_back(e[size() - 1].P); //push the last point in the centerline
|
9c97e126
David Mayerich
added an axis-ali...
|
722
|
|
58ee2b21
David Mayerich
incorporated stim...
|
723
724
|
return cylinder<T>(result);
|
e21d1051
David Mayerich
cylinder can buil...
|
725
|
}*/
|
9c97e126
David Mayerich
added an axis-ali...
|
726
|
|
26aee9ed
Pavel Govyadinov
changed circle cl...
|
727
|
|
5b6c317a
Pavel Govyadinov
implemented the c...
|
728
729
730
|
};
}
|
3b14b0f5
Pavel Govyadinov
added operators f...
|
731
|
#endif
|