File failed to load: https://isis.astrogeology.usgs.gov/9.0.0/Object/assets/jax/output/NativeMML/config.js
Isis 3 Programmer Reference
EmbreeTargetShape.cpp
1
5
6/* SPDX-License-Identifier: CC0-1.0 */
7
8#include "EmbreeTargetShape.h"
9
10#include <iostream>
11#include <iomanip>
12#include <numeric>
13#include <sstream>
14
15#include "NaifDskApi.h"
16
17#include "FileName.h"
18#include "IException.h"
19#include "IString.h"
20#include "NaifStatus.h"
21#include "Pvl.h"
22
23namespace Isis {
24
29 ray.tnear = 0.0;
30 ray.tfar = std::numeric_limits<float>::infinity(); // Should be INF
31 ray.mask = 0xFFFFFFFF;
32 lastHit = -1;
33 hit.u = 0.0;
34 hit.v = 0.0;
35 hit.geomID = RTC_INVALID_GEOMETRY_ID;
36 hit.primID = RTC_INVALID_GEOMETRY_ID;
37 hit.instID[0] = RTC_INVALID_GEOMETRY_ID;
38 ray.org_x = 0.0;
39 ray.org_y = 0.0;
40 ray.org_z = 0.0;
41 ray.dir_x = 0.0;
42 ray.dir_y = 0.0;
43 ray.dir_z = 0.0;
44 hit.Ng_x = 0.0;
45 hit.Ng_y = 0.0;
46 hit.Ng_z = 0.0;
47 }
48
49
56 RTCMultiHitRay::RTCMultiHitRay(const std::vector<double> &origin,
57 const std::vector<double> &direction) {
58 ray.tnear = 0.0;
59 ray.tfar = std::numeric_limits<float>::infinity(); // Should be INF
60 ray.mask = 0xFFFFFFFF;
61 lastHit = -1;
62 hit.u = 0.0;
63 hit.v = 0.0;
64 hit.geomID = RTC_INVALID_GEOMETRY_ID;
65 hit.primID = RTC_INVALID_GEOMETRY_ID;
66 hit.instID[0] = RTC_INVALID_GEOMETRY_ID;
67 ray.org_x = origin[0];
68 ray.org_y = origin[1];
69 ray.org_z = origin[2];
70 ray.dir_x = direction[0];
71 ray.dir_y = direction[1];
72 ray.dir_z = direction[2];
73 hit.Ng_x = 0.0;
74 hit.Ng_y = 0.0;
75 hit.Ng_z = 0.0;
76 }
77
78
87 LinearAlgebra::Vector direction) {
88 ray.tnear = 0.0;
89 ray.tfar = std::numeric_limits<float>::infinity(); // Should be INF
90 ray.mask = 0xFFFFFFFF;
91 lastHit = -1;
92 hit.u = 0.0;
93 hit.v = 0.0;
94 hit.geomID = RTC_INVALID_GEOMETRY_ID;
95 hit.primID = RTC_INVALID_GEOMETRY_ID;
96 hit.instID[0] = RTC_INVALID_GEOMETRY_ID;
97 ray.org_x = origin[0];
98 ray.org_y = origin[1];
99 ray.org_z = origin[2];
100 ray.dir_x = direction[0];
101 ray.dir_y = direction[1];
102 ray.dir_z = direction[2];
103 hit.Ng_x = 0.0;
104 hit.Ng_y = 0.0;
105 hit.Ng_z = 0.0;
106 }
107
108
113 ray.tnear = 0.0;
114 ray.tfar = std::numeric_limits<float>::infinity(); // Should be INF
115 ray.mask = 0xFFFFFFFF;
116 lastHit = -1;
117 hit.u = 0.0;
118 hit.v = 0.0;
119 ignorePrimID = -1;
120 hit.geomID = RTC_INVALID_GEOMETRY_ID;
121 hit.primID = RTC_INVALID_GEOMETRY_ID;
122 hit.instID[0] = RTC_INVALID_GEOMETRY_ID;
123 ray.org_x = 0.0;
124 ray.org_y = 0.0;
125 ray.org_z = 0.0;
126 ray.dir_x = 0.0;
127 ray.dir_y = 0.0;
128 ray.dir_z = 0.0;
129 hit.Ng_x = 0.0;
130 hit.Ng_y = 0.0;
131 hit.Ng_z = 0.0;
132 }
133
134
141 RTCOcclusionRay::RTCOcclusionRay(const std::vector<double> &origin,
142 const std::vector<double> &direction) {
143 ray.tnear = 0.0;
144 ray.tfar = std::numeric_limits<float>::infinity(); // Should be INF
145 ray.mask = 0xFFFFFFFF;
146 lastHit = -1;
147 hit.u = 0.0;
148 hit.v = 0.0;
149 hit.geomID = RTC_INVALID_GEOMETRY_ID;
150 ignorePrimID = -1;
151 hit.instID[0] = RTC_INVALID_GEOMETRY_ID;
152 ray.org_x = origin[0];
153 ray.org_y = origin[1];
154 ray.org_z = origin[2];
155 ray.dir_x = direction[0];
156 ray.dir_y = direction[1];
157 ray.dir_z = direction[2];
158 hit.Ng_x = 0.0;
159 hit.Ng_y = 0.0;
160 hit.Ng_z = 0.0;
161 }
162
171 LinearAlgebra::Vector direction) {
172 ray.tnear = 0.0;
173 ray.tfar = std::numeric_limits<float>::infinity(); // Should be INF
174 ray.mask = 0xFFFFFFFF;
175 lastHit = -1;
176 hit.u = 0.0;
177 hit.v = 0.0;
178 hit.geomID = RTC_INVALID_GEOMETRY_ID;
179 ignorePrimID = -1;
180 hit.instID[0] = RTC_INVALID_GEOMETRY_ID;
181 ray.org_x = origin[0];
182 ray.org_y = origin[1];
183 ray.org_z = origin[2];
184 ray.dir_x = direction[0];
185 ray.dir_y = direction[1];
186 ray.dir_z = direction[2];
187 hit.Ng_x = 0.0;
188 hit.Ng_y = 0.0;
189 hit.Ng_z = 0.0;
190 }
191
192
201
202
212 LinearAlgebra::Vector &normal,
213 int prim)
214 : intersection(location),
215 surfaceNormal(normal) ,
216 primID(prim) { }
217
218
226 : m_name(),
227 m_mesh(),
228 m_cloud(),
229 m_device(rtcNewDevice(NULL)),
230 m_scene(rtcNewScene(m_device))
231 {
232 rtcSetSceneFlags(m_scene, RTC_SCENE_FLAG_ROBUST | RTC_SCENE_FLAG_CONTEXT_FILTER_FUNCTION);
233 rtcSetSceneBuildQuality(m_scene, RTC_BUILD_QUALITY_HIGH);
234 }
235
242 EmbreeTargetShape::EmbreeTargetShape(pcl::PolygonMesh::Ptr mesh, const QString &name)
243 : m_name(name),
244 m_mesh(),
245 m_cloud(),
246 m_device(rtcNewDevice(NULL)),
247 m_scene(rtcNewScene(m_device))
248 {
249 rtcSetSceneFlags(m_scene, RTC_SCENE_FLAG_ROBUST | RTC_SCENE_FLAG_CONTEXT_FILTER_FUNCTION);
250 rtcSetSceneBuildQuality(m_scene, RTC_BUILD_QUALITY_HIGH);
251 initMesh(mesh);
252 }
253
264 EmbreeTargetShape::EmbreeTargetShape(const QString &dem, const Pvl *conf)
265 : m_name(),
266 m_mesh(),
267 m_cloud(),
268 m_device(rtcNewDevice(NULL)),
269 m_scene(rtcNewScene(m_device)) {
270 rtcSetSceneFlags(m_scene, RTC_SCENE_FLAG_ROBUST | RTC_SCENE_FLAG_CONTEXT_FILTER_FUNCTION);
271 rtcSetSceneBuildQuality(m_scene, RTC_BUILD_QUALITY_HIGH);
272 FileName file(dem);
273 pcl::PolygonMesh::Ptr mesh;
274 m_name = file.baseName();
275
276 try {
277 // DEMs (ISIS cubes) TODO implement this
278 if (file.extension() == "cub") {
279 QString msg = "DEMs cannot be used to create an EmbreeTargetShape.";
280 throw IException(IException::Io, msg, _FILEINFO_);
281 }
282 // DSKs
283 else if (file.extension().toLower() == "bds") {
284 mesh = readDSK(file);
285 }
286 // Let PCL try to handle other formats (obj, ply, etc.)
287 else {
288 mesh = readPC(file);
289 }
290 }
291 catch (IException &e) {
292 QString msg = "Failed creating an EmbreeTargetShape from ["
293 + file.expanded() + "].";
294 throw IException(e, IException::Io, msg, _FILEINFO_);
295 }
296 initMesh(mesh);
297 }
298
299
311 pcl::PolygonMesh::Ptr EmbreeTargetShape::readDSK(FileName file) {
312
314 SpiceInt dskHandle;
315 SpiceDLADescr dlaDescriptor;
317 SpiceDSKDescr dskDescriptor;
318 SpiceInt numPlates;
319 SpiceInt numVertices;
320
321 // Sanity check
322 if ( !file.fileExists() ) {
323 QString mess = "NAIF DSK file [" + file.expanded() + "] does not exist.";
324 throw IException(IException::User, mess, _FILEINFO_);
325 }
326
327 // Open the NAIF Digital Shape Kernel (DSK)
328 dasopr_c( file.expanded().toLatin1().data(), &dskHandle );
330
331 // Search to the first DLA segment
332 SpiceBoolean found;
333 dlabfs_c( dskHandle, &dlaDescriptor, &found );
335 if ( !found ) {
336 QString mess = "No segments found in DSK file [" + file.expanded() + "]";
337 throw IException(IException::User, mess, _FILEINFO_);
338 }
339
340 dskgd_c( dskHandle, &dlaDescriptor, &dskDescriptor );
342
343 // Get The number of polygons and vertices
344 dskz02_c( dskHandle, &dlaDescriptor, &numVertices, &numPlates );
346
347 // Allocate polygon and vertices arrays
348 // These arrays are actually numVertices x 3 and numPlates x 3,
349 // this allocation is easier and produces needed memory layout for dskv02_c
350 // These arrays are very large, so they need to be heap allocated.
351 double *verticesArray = new double[numVertices * 3];
352 int *polygonsArray = new int[numPlates * 3];
353
354 // Read the vertices from the dsk file
355 SpiceInt numRead = 0;
356 dskv02_c(dskHandle, &dlaDescriptor, 1, numVertices,
357 &numRead, ( SpiceDouble(*)[3] )(verticesArray) );
359 if ( numRead != numVertices ) {
360 QString msg = "Failed reading all vertices from the DSK file, ["
361 + toString(numRead) + "] out of ["
362 + toString(numVertices) + "] vertices read.";
363 throw IException(IException::Io, msg, _FILEINFO_);
364 }
365
366 // Read the polygons from the DSK
367 numRead = 0;
368 dskp02_c(dskHandle, &dlaDescriptor, 1, numPlates,
369 &numRead, ( SpiceInt(*)[3] )(polygonsArray) );
371 if ( numRead != numPlates ) {
372 QString msg = "Failed reading all polygons from the DSK file, ["
373 + toString(numRead) + "] out of ["
374 + toString(numPlates) + "] polygons read.";
375 throw IException(IException::Io, msg, _FILEINFO_);
376 }
377
378 // Store the vertices in a point cloud
379 pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
380 for (int vertexIndex = 0; vertexIndex < numVertices; ++vertexIndex) {
381 cloud->push_back(pcl::PointXYZ(verticesArray[vertexIndex * 3],
382 verticesArray[vertexIndex * 3 + 1],
383 verticesArray[vertexIndex * 3 + 2]));
384 }
385
386 // Store the polygons as a vector of vertices
387 std::vector<pcl::Vertices> polygons;
388 for (int plateIndex = 0; plateIndex < numPlates; ++plateIndex) {
389 pcl::Vertices vertexIndices;
390 // NAIF uses 1 based indexing for the vertices, so subtract 1
391 vertexIndices.vertices.push_back(polygonsArray[plateIndex * 3] - 1);
392 vertexIndices.vertices.push_back(polygonsArray[plateIndex * 3 + 1] - 1);
393 vertexIndices.vertices.push_back(polygonsArray[plateIndex * 3 + 2] - 1);
394 polygons.push_back(vertexIndices);
395 }
396
397 // Create the mesh
398 pcl::PolygonMesh::Ptr mesh(new pcl::PolygonMesh);
399 mesh->polygons = polygons;
400 pcl::PCLPointCloud2::Ptr cloudBlob(new pcl::PCLPointCloud2);
401 pcl::toPCLPointCloud2(*cloud, *cloudBlob);
402 mesh->cloud = *cloudBlob;
403
404 // Free the vectors used to read the vertices and polygons
405 delete [] verticesArray;
406 delete [] polygonsArray;
407
408 return mesh;
409 }
410
411
421 pcl::PolygonMesh::Ptr EmbreeTargetShape::readPC(FileName file) {
422 pcl::PolygonMesh::Ptr mesh( new pcl::PolygonMesh );
423
424 int loadStatus = pcl::io::load(file.expanded().toStdString(), *mesh);
425 if (loadStatus == -1) {
426 QString msg = "Failed loading target shape file [" + file.expanded() + "]";
427 throw IException(IException::Io, msg, _FILEINFO_);
428 }
429 return mesh;
430 }
431
432
450 void EmbreeTargetShape::initMesh(pcl::PolygonMesh::Ptr mesh) {
451 m_mesh = mesh;
452
453 // The points are stored in a pcl::PCLPointCloud2 object that we cannot used.
454 // So, convert them into a pcl::PointCloud<pcl::PointXYZ> object that we can use.
455 pcl::fromPCLPointCloud2(mesh->cloud, m_cloud);
456
457 // Create a static geometry (the body) in our scene
458 RTCGeometry rtcMesh = rtcNewGeometry(m_device, RTC_GEOMETRY_TYPE_TRIANGLE);
459 // Add vertices and faces (indices)
460 if (!isValid()) {
461 return;
462 }
463 // Add the body's vertices to the Embree ray tracing device's vertex buffer
464 Vertex *vertices = (Vertex *)rtcSetNewGeometryBuffer(rtcMesh, RTC_BUFFER_TYPE_VERTEX, 0, RTC_FORMAT_FLOAT3, sizeof(Vertex), numberOfVertices());
465 for (int v = 0; v < numberOfVertices(); ++v) {
466 vertices[v].x = m_cloud.points[v].x;
467 vertices[v].y = m_cloud.points[v].y;
468 vertices[v].z = m_cloud.points[v].z;
469 }
470
471 if (!isValid()) {
472 return;
473 }
474 Triangle *triangles = (Triangle *)rtcSetNewGeometryBuffer(rtcMesh, RTC_BUFFER_TYPE_INDEX, 0, RTC_FORMAT_UINT3, sizeof(Triangle), numberOfPolygons());
475 // Add the body's face (vertex indices) to the Embree device's index buffer
476 for (int t = 0; t < numberOfPolygons(); ++t) {
477 triangles[t].v0 = m_mesh->polygons[t].vertices[0];
478 triangles[t].v1 = m_mesh->polygons[t].vertices[1];
479 triangles[t].v2 = m_mesh->polygons[t].vertices[2];
480 }
481
482 // Add the multi-hit filter
483 rtcSetGeometryIntersectFilterFunction(rtcMesh, EmbreeTargetShape::multiHitFilter);
484
485 // Add the occlusion filter
486 rtcSetGeometryOccludedFilterFunction(rtcMesh, EmbreeTargetShape::occlusionFilter);
487
488 rtcCommitGeometry(rtcMesh);
489 rtcAttachGeometry(m_scene, rtcMesh);
490 rtcReleaseGeometry(rtcMesh);
491
492 // Done, now we can perform some ray tracing
493 rtcCommitScene(m_scene);
494 }
495
496
507 // if (!isValid()) {
508 // return;
509 // }
510 // // Add the body's vertices to the Embree ray tracing device's vertex buffer
511 // Vertex *vertices = (Vertex *) rtcMapBuffer(m_scene, geomID, RTC_VERTEX_BUFFER);
512 // for (int v = 0; v < numberOfVertices(); ++v) {
513 // vertices[v].x = m_cloud.points[v].x;
514 // vertices[v].y = m_cloud.points[v].y;
515 // vertices[v].z = m_cloud.points[v].z;
516 // }
517 // // Flush buffer
518 // rtcUnmapBuffer(m_scene, geomID, RTC_VERTEX_BUFFER);
519 }
520
521
532 // if (!isValid()) {
533 // return;
534 // }
535 // // Add the body's face (vertex indices) to the Embree device's index buffer
536 // Triangle *triangles = (Triangle *) rtcMapBuffer(m_scene, geomID, RTC_INDEX_BUFFER);
537 // for (int t = 0; t < numberOfPolygons(); ++t) {
538 // triangles[t].v0 = m_mesh->polygons[t].vertices[0];
539 // triangles[t].v1 = m_mesh->polygons[t].vertices[1];
540 // triangles[t].v2 = m_mesh->polygons[t].vertices[2];
541 // }
542 // // Flush buffer
543 // rtcUnmapBuffer(m_scene, geomID, RTC_INDEX_BUFFER);
544 }
545
546
552 rtcReleaseScene(m_scene);
553 rtcReleaseDevice(m_device);
554 }
555
556
564 if (isValid()) {
565 return m_mesh->polygons.size();
566 }
567 return 0;
568 }
569
570
578 if (isValid()) {
579 return m_mesh->cloud.height * m_mesh->cloud.width;
580 }
581 return 0;
582 }
583
584
596 RTCBounds sceneBounds;
597 if (isValid()) {
598 rtcGetSceneBounds(m_scene, &sceneBounds);
599 }
600 else {
601 sceneBounds.lower_x = 0.0;
602 sceneBounds.lower_y = 0.0;
603 sceneBounds.lower_z = 0.0;
604 sceneBounds.upper_x = 0.0;
605 sceneBounds.upper_y = 0.0;
606 sceneBounds.upper_z = 0.0;
607 }
608 return sceneBounds;
609 }
610
611
619 RTCBounds sceneBoundary = sceneBounds();
620 LinearAlgebra::Vector diagonal(3);
621 diagonal[0] = sceneBoundary.upper_x - sceneBoundary.lower_x;
622 diagonal[1] = sceneBoundary.upper_y - sceneBoundary.lower_y;
623 diagonal[2] = sceneBoundary.upper_z - sceneBoundary.lower_z;
624 return LinearAlgebra::magnitude(diagonal);
625 }
626
627
658 if (isValid()) {
659 RTCIntersectContext context;
660 rtcInitIntersectContext(&context);
661
662 rtcIntersect1(m_scene, &context, (RTCRayHit *)&ray);
663 }
664 }
665
666
675 RTCIntersectContext context;
676 rtcInitIntersectContext(&context);
677
678 rtcOccluded1(m_scene, &context, (RTCRay*)&ray);
679
680 // rtcOccluded sets the ray.tfar to -inf if the ray hits anything
681 if (isinf(ray.ray.tfar) && ray.ray.tfar < 0) {
682 return true;
683 }
684 return false;
685 }
686
687
708 if (hitIndex > ray.lastHit || hitIndex < 0) {
709 QString msg = "Hit index [" + toString(hitIndex) + "] is out of bounds.";
710 throw IException(IException::Programmer, msg, _FILEINFO_);
711 }
712 RTCGeometry geom = rtcGetGeometry(m_scene, ray.hitGeomIDs[hitIndex]);
713 Vertex *vertices = (Vertex *)rtcGetGeometryBufferData(geom, RTC_BUFFER_TYPE_VERTEX, 0);
714 Triangle *triangles = (Triangle *)rtcGetGeometryBufferData(geom, RTC_BUFFER_TYPE_INDEX, 0);
715
716 // Get the vertices of the triangle hit
717 Vertex v0 = vertices[triangles[ray.hitPrimIDs[hitIndex]].v0];
718 Vertex v1 = vertices[triangles[ray.hitPrimIDs[hitIndex]].v1];
719 Vertex v2 = vertices[triangles[ray.hitPrimIDs[hitIndex]].v2];
720
721 // The intersection location comes out in barycentric coordinates, (u, v, w).
722 // Only u and v are returned because u + v + w = 1. If the coordinates of the
723 // triangle vertices are v0, v1, and v2, then the cartesian coordinates are:
724 // w*v0 + u*v1 + v*v2
725 float u = ray.hitUs[hitIndex];
726 float v = ray.hitVs[hitIndex];
727 float w = 1.0 - u - v;
728
729 LinearAlgebra::Vector intersection(3);
730 intersection[0] = w*v0.x + u*v1.x + v*v2.x;
731 intersection[1] = w*v0.y + u*v1.y + v*v2.y;
732 intersection[2] = w*v0.z + u*v1.z + v*v2.z;
733
734 // Calculate the normal vector as (v1 - v0) x (v2 - v0) and normalize it
735 // TODO This calculation assumes that the shape conforms to the NAIF dsk standard
736 // of the plate vertices being ordered counterclockwise about the normal.
737 // Check if this is true for other file types and/or make a more generic process
738 LinearAlgebra::Vector surfaceNormal(3);
739 surfaceNormal[0] = (v1.y - v0.y) * (v2.z - v0.z)
740 - (v1.z - v0.z) * (v2.y - v0.y);
741 surfaceNormal[1] = (v1.z - v0.z) * (v2.x - v0.x)
742 - (v1.x - v0.x) * (v2.z - v0.z);
743 surfaceNormal[2] = (v1.x - v0.x) * (v2.y - v0.y)
744 - (v1.y - v0.y) * (v2.x - v0.x);
745
746 // The surface normal is not normalized so normalize it.
747 surfaceNormal = LinearAlgebra::normalize(surfaceNormal);
748 return RayHitInformation(intersection, surfaceNormal, ray.hitPrimIDs[hitIndex]);
749 }
750
751
757 QString EmbreeTargetShape::name() const {
758 return ( m_name );
759 }
760
761
768 return m_mesh.get();
769 }
770
771
781 void EmbreeTargetShape::multiHitFilter(const RTCFilterFunctionNArguments* args) {
782 /* avoid crashing when debug visualizations are used */
783 if (args->context == nullptr)
784 return;
785
786 assert(args->N == 1);
787 int *valid = args->valid;
788 if (valid[0] != -1) {
789 return;
790 }
791 RTCMultiHitRay *ray = (RTCMultiHitRay *)args->ray;
792 RTCHit *hit = (RTCHit *)args->hit;
793
794 // Calculate the index to store the hit in
795 ray->lastHit++;
796
797 // Store the hits
798 ray->hitGeomIDs[ray->lastHit] = hit->geomID;
799 ray->hitPrimIDs[ray->lastHit] = hit->primID;
800 ray->hitUs[ray->lastHit] = hit->u;
801 ray->hitVs[ray->lastHit] = hit->v;
802
803 // If there are less than 16 hits, continue ray tracing.
804 if (ray->lastHit < 15) {
805 valid[0] = 0;
806 ray->hit.geomID = RTC_INVALID_GEOMETRY_ID;
807 hit->geomID = RTC_INVALID_GEOMETRY_ID;
808 }
809 }
810
820 void EmbreeTargetShape::occlusionFilter(const RTCFilterFunctionNArguments* args) {
821 /* avoid crashing when debug visualizations are used */
822 if (args->context == nullptr)
823 return;
824
825 assert(args->N == 1);
826 int *valid = args->valid;
827 if (valid[0] != -1) {
828 return;
829 }
830 RTCOcclusionRay *ray = (RTCOcclusionRay *)args->ray;
831 RTCHit *hit = (RTCHit *)args->hit;
832 ray->hit.primID = hit->primID;
833
834 // This is the case where we've re-intersected the occluded plate. If this happens, ignore
835 // and keep tracing
836 if (ray->hit.primID == ray->ignorePrimID) {
837 valid[0] = 0;
838 // ray->hit.geomID = RTC_INVALID_GEOMETRY_ID;
839 // hit->geomID = RTC_INVALID_GEOMETRY_ID;
840 }
841 }
842
843} // namespace Isis
bool isOccluded(RTCOcclusionRay &ray)
Check if a ray intersects the target body.
void initMesh(pcl::PolygonMesh::Ptr mesh)
Internalize a PointCloudLibrary polygon mesh in the target shape.
pcl::PolygonMesh::Ptr readDSK(FileName file)
Read a NAIF type 2 DSK file into a PointCloudLibrary polygon mesh.
static void occlusionFilter(const RTCFilterFunctionNArguments *args)
Filter function for collecting multiple primitiveIDs This function is called by the Embree library du...
int numberOfPolygons() const
Return the number of polygons in the target shape.
int numberOfVertices() const
Return the number of vertices in the target shape.
virtual ~EmbreeTargetShape()
Desctructor.
bool isValid() const
Return if a valid mesh is internalized and ready for use.
RTCDevice m_device
!< The point cloud representation of the target.
void addVertices(int geomID)
Adds the vertices from the internalized vertex point cloud to the Embree scene.
static void multiHitFilter(const RTCFilterFunctionNArguments *args)
Filter function for collecting multiple hits during ray intersection.
pcl::PolygonMesh::Ptr readPC(FileName file)
Read a PointCloudLibrary file into a PointCloudLibrary polygon mesh.
pcl::PolygonMesh::Ptr m_mesh
!< The name of the target.
RTCScene m_scene
!< The Embree device for rendering the scene.
void addIndices(int geomID)
Adds the polygon vertex indices from the internalized polygon mesh to the Embree scene.
RayHitInformation getHitInformation(RTCMultiHitRay &ray, int hitIndex)
Extract the intersection point and unit surface normal from an RTCMultiHitRay that has been intersect...
RTCBounds sceneBounds() const
Returns the bounds of the Embree scene.
void intersectRay(RTCMultiHitRay &ray)
Intersect a ray with the target shape.
EmbreeTargetShape()
Default empty constructor.
pcl::PointCloud< pcl::PointXYZ > m_cloud
!< A boost shared pointer to the polygon mesh representation of the target.
QString name() const
Return the name of the target shape.
double maximumSceneDistance() const
Return the maximum distance within the scene.
File name manipulation and expansion.
Definition FileName.h:100
bool fileExists() const
Returns true if the file exists; false otherwise.
Definition FileName.cpp:449
QString baseName() const
Returns the name of the file without the path and without extensions.
Definition FileName.cpp:145
QString expanded() const
Returns a QString of the full file name including the file path, excluding the attributes.
Definition FileName.cpp:196
QString extension() const
Returns the last extension of the file name.
Definition FileName.cpp:178
Isis exception class.
Definition IException.h:91
@ User
A type of error that could only have occurred due to a mistake on the user's part (e....
Definition IException.h:126
@ Programmer
This error is for when a programmer made an API call that was illegal.
Definition IException.h:146
@ Io
A type of error that occurred when performing an actual I/O operation.
Definition IException.h:155
static Vector vector(double v0, double v1, double v2)
Constructs a 3 dimensional vector with the given component values.
static double magnitude(const Vector &vector)
Computes the magnitude (i.e., the length) of the given vector using the Euclidean norm (L2 norm).
boost::numeric::ublas::vector< double > Vector
Definition for an Isis::LinearAlgebra::Vector of doubles.
static Vector normalize(const Vector &vector)
Returns a unit vector that is codirectional with the given vector by dividing each component of the v...
static void CheckErrors(bool resetNaif=true)
This method looks for any naif errors that might have occurred.
Container for cube-like labels.
Definition Pvl.h:119
This is free and unencumbered software released into the public domain.
Definition Apollo.h:16
QString toString(bool boolToConvert)
Global function to convert a boolean to a string.
Definition IString.cpp:211
Container for a tin, or triangular polygon.
int v1
The index of the second vertex in the tin.
int v2
The index of the third vertex in the tin.
int v0
The index of the first vertex in the tin.
Struct for capturing multiple intersections when using embree::rtcintersectscene.
float hitUs[16]
Barycentric u coordinate of the hits.
unsigned hitGeomIDs[16]
IDs of the geometries (bodies) hit.
float hitVs[16]
Barycentric v coordinate of the hits.
unsigned hitPrimIDs[16]
IDs of the primitives (trinagles) hit.
RTCMultiHitRay()
Default constructor for RTCMultiHitRay.
int lastHit
Index of the last hit in the hit containers.
Struct for capturing occluded plates when using embree::rtcintersectscene.
RTCOcclusionRay()
Default constructor for RTCOcclussionRay.
unsigned ignorePrimID
IDs of the primitives (trinagles) which should be ignored.
int lastHit
Index of the last hit in the hit containers.
Container that holds the body fixed intersection point and unit surface normal for a hit.
int primID
The primitive ID of the hit.
LinearAlgebra::Vector intersection
The (x, y, z) intersection location.
LinearAlgebra::Vector surfaceNormal
The unit surface normal vector at the intersection.
RayHitInformation()
Default constructor for RayHitInformation.