23 #ifndef STLFILEWRITER_HH 24 #define STLFILEWRITER_HH 32 #include <blitz/array.h> 39 std::ofstream _outFile;
40 std::string _objectName;
45 int _nEllipseSegments;
50 const std::string& objectName=
"dummy")
51 : _outFile(fileName.c_str()),
52 _objectName(objectName),
59 _outFile <<
"solid " << _objectName <<
"\n";
64 _outFile <<
"endsolid " << _objectName <<
"\n";
70 const blitz::TinyVector<float,3>& p2,
71 const blitz::TinyVector<float,3>& p3,
72 const blitz::TinyVector<float,3>& normal =
73 blitz::TinyVector<float,3>(0,0,0))
75 _outFile <<
"facet normal " 76 << normal(0) <<
" " << normal(1) <<
" " << normal(2)
80 << p1(0) <<
" " << p1(1) <<
" " << p1(2) <<
"\n" 82 << p2(0) <<
" " << p2(1) <<
" " << p2(2) <<
"\n" 84 << p3(0) <<
" " << p3(1) <<
" " << p3(2) <<
"\n" 90 const blitz::TinyVector<float,3>& origin,
91 const blitz::TinyVector<float,3>& nx,
92 const blitz::TinyVector<float,3>& ny,
93 const blitz::TinyVector<float,3>& nz,
94 const std::vector< blitz::TinyVector<float,2> >& contour,
95 const std::vector< blitz::TinyVector<float,2> >& normals,
99 for(
int i = 0; i < nAngles; ++i)
101 float phi1 = float(i)/nAngles*2*M_PI;
102 float phi2 = float(i+1)/nAngles*2*M_PI;
103 blitz::TinyVector<float,3> rotNy1 =
104 sin( phi1) * ny + cos( phi1) * nz;
106 blitz::TinyVector<float,3> rotNy2 =
107 sin( phi2) * ny + cos( phi2) * nz;
109 blitz::TinyVector<float,3> rotNyMean =
111 rotNyMean /= sqrt(
blitz::dot( rotNyMean, rotNyMean));
113 for(
size_t j = 0; j < contour.size()-1; ++j)
115 blitz::TinyVector<float,3> p1 =
118 + rotNy1 * contour[j](1);
120 blitz::TinyVector<float,3> p2 =
123 + rotNy2 * contour[j](1);
125 blitz::TinyVector<float,3> p3 =
127 + nx * contour[j+1](0)
128 + rotNy1 * contour[j+1](1);
130 blitz::TinyVector<float,3> p4 =
132 + nx * contour[j+1](0)
133 + rotNy2 * contour[j+1](1);
135 blitz::TinyVector<float,3> normal =
137 + rotNyMean * normals[j](1);
157 const blitz::TinyVector<float,3>& dxyz,
160 if (radius <= 0) radius = _tailRadius;
163 float length = sqrt(
blitz::dot( dxyz, dxyz));
165 blitz::TinyVector<float,3> nx = dxyz / length;
167 blitz::TinyVector<float,3> unitCol(0,0,1);
168 blitz::TinyVector<float,3> ny;
169 blitz::TinyVector<float,3> nz;
171 ny = blitz::cross( nx, unitCol);
172 nz = blitz::cross( nx, ny);
178 std::vector<blitz::TinyVector<float,2> > contour(4);
179 contour[0] = length, 0;
180 contour[1] = length, radius;
181 contour[2] = 0, radius;
184 std::vector<blitz::TinyVector<float,2> > normals(4);
190 for(
size_t i = 0; i < normals.size(); ++i)
192 normals[i] /= sqrt(
blitz::dot( normals[i], normals[i]));
197 contour, normals, _nSegments);
202 void drawArrow(
const blitz::TinyVector<float,3>& origin,
203 const blitz::TinyVector<float,3>& dxyz)
206 float length = std::sqrt(
blitz::dot( dxyz, dxyz));
208 blitz::TinyVector<float,3> nx = dxyz / length;
210 blitz::TinyVector<float,3> unitCol(0,0,1);
211 blitz::TinyVector<float,3> ny;
212 blitz::TinyVector<float,3> nz;
214 ny = blitz::cross( nx, unitCol);
215 nz = blitz::cross( nx, ny);
221 std::vector< blitz::TinyVector<float,2> > contour(5);
222 contour[0] = length, 0;
223 contour[1] = length-2*_headRadius, _headRadius;
224 contour[2] = length-2*_headRadius, _tailRadius;
225 contour[3] = 0, _tailRadius;
228 std::vector< blitz::TinyVector<float,2> > normals(5);
235 for(
size_t i = 0; i < normals.size(); ++i)
237 normals[i] /= std::sqrt(
blitz::dot( normals[i], normals[i]));
242 contour, normals, _nSegments);
247 const blitz::TinyVector<float,3>& axis1,
248 const blitz::TinyVector<float,3>& axis2)
250 float radius1 = sqrt(
blitz::dot( axis1, axis1));
251 float radius2 = sqrt(
blitz::dot( axis2, axis2));
253 blitz::TinyVector<float,3> ny = axis1 / radius1;
254 blitz::TinyVector<float,3> nz = axis2 / radius2;
255 blitz::TinyVector<float,3> nx;
256 nx = blitz::cross( ny, nz);
259 std::vector<blitz::TinyVector<float,2> > contour(_nSegments);
260 std::vector<blitz::TinyVector<float,2> > normals(_nSegments);
261 for(
int i = 0; i < _nSegments; ++i)
263 float phi = float(i)/(_nSegments-1)*2*M_PI;
265 _borderRadius*sin(phi),
266 _borderRadius*cos(phi);
268 normals[i] = sin(phi), cos(phi);
271 for(
size_t i = 0; i < normals.size(); ++i)
273 normals[i] /= sqrt(
blitz::dot( normals[i], normals[i]));
277 for(
int i = 0; i < _nEllipseSegments; ++i)
279 float phi1 = float(i)/_nEllipseSegments*2*M_PI;
280 float phi2 = float(i+1)/_nEllipseSegments*2*M_PI;
281 blitz::TinyVector<float,3> center1 =
282 radius1 * sin( phi1) * ny + radius2 * cos( phi1) * nz;
284 blitz::TinyVector<float,3> rotNy1 =
285 center1 / sqrt(
blitz::dot( center1, center1));
287 blitz::TinyVector<float,3> center2 =
288 radius1 * sin( phi2) * ny + radius2 * cos( phi2) * nz;
290 blitz::TinyVector<float,3> rotNy2 =
291 center2 / sqrt(
blitz::dot( center2, center2));
293 blitz::TinyVector<float,3> rotNyMean =
295 rotNyMean /= sqrt(
blitz::dot( rotNyMean, rotNyMean));
297 for(
size_t j = 0; j < contour.size()-1; ++j)
299 blitz::TinyVector<float,3> p1 =
302 + center1 + rotNy1 * contour[j](1);
304 blitz::TinyVector<float,3> p2 =
307 + center2 + rotNy2 * contour[j](1);
309 blitz::TinyVector<float,3> p3 =
311 + nx * contour[j+1](0)
312 + center1 + rotNy1 * contour[j+1](1);
314 blitz::TinyVector<float,3> p4 =
316 + nx * contour[j+1](0)
317 + center2 + rotNy2 * contour[j+1](1);
319 blitz::TinyVector<float,3> normal =
321 + rotNyMean * normals[j](1);
342 const float radius,
const int halfCircleSampling = 10)
345 const int sampling = halfCircleSampling;
346 blitz::TinyVector<float,3> znormal(0,0,1);
348 blitz::Array< blitz::TinyVector<float,3>, 1>
349 startHalfCircleSamples(sampling);
350 startHalfCircleSamples = blitz::TinyVector<float,3>(0,0,0);
352 for (
int i = 0; i < sampling; i++)
354 float factor = ( float(i)/float(sampling-1)) - 0.5;
355 float x = radius * cos( factor * M_PI);
357 float z = radius * sin( factor * M_PI);
358 startHalfCircleSamples(i) = blitz::TinyVector<float,3>(x,y,z);
361 blitz::Array< blitz::TinyVector<float,3>, 1>
362 halfCircleSamples(sampling);
363 halfCircleSamples = startHalfCircleSamples;
364 blitz::Array< blitz::TinyVector<float,3>, 1>
365 nextHalfCircleSamples(sampling);
366 nextHalfCircleSamples = blitz::TinyVector<float,3>(0,0,0);
368 for (
int i = 1; i < 2*sampling; ++i)
370 float factor = ( float(i)/float(2*sampling-1));
371 float phi = factor * 2*M_PI;
372 float cosPhi = cos(phi);
373 float sinPhi = sin(phi);
376 blitz::Array<float,2> rotMatrix(3,3);
377 rotMatrix(0,0) = cosPhi;
378 rotMatrix(0,1) = -sinPhi;
381 rotMatrix(1,0) = sinPhi;
382 rotMatrix(1,1) = cosPhi;
389 for(
int j = 0; j < startHalfCircleSamples.extent(0); ++j) {
390 nextHalfCircleSamples(j) =
391 mvMult(rotMatrix, startHalfCircleSamples(j));
394 for (
int k = 0; k < sampling-2; ++k)
396 blitz::TinyVector<float,3> p1 = nextHalfCircleSamples(k);
397 blitz::TinyVector<float,3> p2 = nextHalfCircleSamples(k+1);
398 blitz::TinyVector<float,3> p3 = halfCircleSamples(k+1);
399 blitz::TinyVector<float,3> v1 = p2 - p1;
400 blitz::TinyVector<float,3> v2 = p3 - p1;
401 blitz::TinyVector<float,3> n = blitz::cross(v1,v2);
402 n /= sqrt(n(0)*n(0) + n(1)*n(1) + n(2)*n(2));
403 p1 += center; p2 += center; p3 += center;
407 for (
int l = 1; l < sampling-1; ++l)
409 blitz::TinyVector<float,3> p1 = halfCircleSamples(l);
410 blitz::TinyVector<float,3> p2 = nextHalfCircleSamples(l);
411 blitz::TinyVector<float,3> p3 = halfCircleSamples(l+1);
412 blitz::TinyVector<float,3> v1 = p2 - p1;
413 blitz::TinyVector<float,3> v2 = p3 - p1;
414 blitz::TinyVector<float,3> n = blitz::cross(v1,v2);
415 n /= sqrt(n(0)*n(0) + n(1)*n(1) + n(2)*n(2));
416 p1 += center; p2 += center; p3 += center;
420 halfCircleSamples = nextHalfCircleSamples;
427 blitz::TinyVector<double,3>
const ¢er,
428 blitz::Array<double,2>
const &f)
430 for (ptrdiff_t thetaIdx = 0; thetaIdx < f.extent(0) - 1; ++thetaIdx)
433 static_cast<double>(thetaIdx) /
434 static_cast<double>(f.extent(0) - 1) * M_PI;
436 static_cast<double>(thetaIdx + 1) /
437 static_cast<double>(f.extent(0) - 1) * M_PI;
438 for (ptrdiff_t phiIdx = 0; phiIdx < f.extent(1); ++phiIdx)
440 double phi =
static_cast<double>(phiIdx) /
441 static_cast<double>(f.extent(1)) * 2.0 * M_PI;
442 double phip = (phiIdx == f.extent(1) - 1) ?
443 0.0 : (static_cast<double>(phiIdx + 1) /
444 static_cast<double>(f.extent(1)) * 2.0 * M_PI);
446 double r = f(thetaIdx, phiIdx);
447 blitz::TinyVector<double,3> v1(
449 r * std::sin(theta) * std::sin(phi),
450 r * std::sin(theta) * std::cos(phi));
454 (phiIdx == f.extent(1) - 1) ? 0 : (phiIdx + 1));
455 blitz::TinyVector<double,3> v2(
456 r * std::cos(thetap),
457 r * std::sin(thetap) * std::sin(phip),
458 r * std::sin(thetap) * std::cos(phip));
461 r = f(thetaIdx + 1, phiIdx);
462 blitz::TinyVector<double,3> v3(
463 r * std::cos(thetap),
464 r * std::sin(thetap) * std::sin(phi),
465 r * std::sin(thetap) * std::cos(phi));
469 (phiIdx == f.extent(1) - 1) ? 0 : (phiIdx + 1));
470 blitz::TinyVector<double,3> v4(
472 r * std::sin(theta) * std::sin(phip),
473 r * std::sin(theta) * std::cos(phip));
476 blitz::TinyVector<double,3> d21(v2 - v1);
477 blitz::TinyVector<double,3> d31(v3 - v1);
478 blitz::TinyVector<double,3> d41(v4 - v1);
480 blitz::TinyVector<double,3> n(
void drawArrow(const blitz::TinyVector< float, 3 > &origin, const blitz::TinyVector< float, 3 > &dxyz)
void drawTriangle(const blitz::TinyVector< float, 3 > &p1, const blitz::TinyVector< float, 3 > &p2, const blitz::TinyVector< float, 3 > &p3, const blitz::TinyVector< float, 3 > &normal=blitz::TinyVector< float, 3 >(0, 0, 0))
void drawEllipse(const blitz::TinyVector< float, 3 > &origin, const blitz::TinyVector< float, 3 > &axis1, const blitz::TinyVector< float, 3 > &axis2)
bool any(blitz::TinyMatrix< bool, NRows, NColumns > const &matrix)
any() reduction for boolean blitz::TinyMatrix.
void drawSphere(const blitz::TinyVector< float, 3 > ¢er, const float radius, const int halfCircleSampling=10)
void drawRotationObjectFromContour(const blitz::TinyVector< float, 3 > &origin, const blitz::TinyVector< float, 3 > &nx, const blitz::TinyVector< float, 3 > &ny, const blitz::TinyVector< float, 3 > &nz, const std::vector< blitz::TinyVector< float, 2 > > &contour, const std::vector< blitz::TinyVector< float, 2 > > &normals, const int nAngles)
blitz::TinyVector< BaseT, Dim > mvMult(MatrixT const &m, blitz::TinyVector< BaseT, Dim > const &v)
Quadratic Matrix * Vector product (generic).
void drawStarshapedSurface(blitz::TinyVector< double, 3 > const ¢er, blitz::Array< double, 2 > const &f)
STLFileWriter(const std::string &fileName, const std::string &objectName="dummy")
void drawCylinder(const blitz::TinyVector< float, 3 > &origin, const blitz::TinyVector< float, 3 > &dxyz, double radius=-1)
double dot(blitz::Array< DataT, 1 > const &vec1, blitz::Array< DataT, 1 > const &vec2)
Calculation of the inner product (dot product) of the given vectors.