Mercurial > hg > orthanc-stone

comparison OrthancStone/Sources/Toolbox/DicomStructureSet.cpp @ 1904:738814c24574

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
speed up rendering of axial slices of RT-STRUCT
author Sebastien Jodogne <s.jodogne@gmail.com>
date Mon, 31 Jan 2022 18:22:57 +0100
parents a5e54bd87b25
children e318b524ad3f
comparison
equal deleted inserted replaced
1903:c7bc9e1776a6 1904:738814c24574
272 return true; 272 return true;
273 } 273 }
274 } 274 }
275 } 275 }
276 276
277 bool DicomStructureSet::Polygon::IsOnSlice(const CoordinateSystem3D& slice) const 277 bool DicomStructureSet::Polygon::IsOnSlice(const CoordinateSystem3D& slice,
278 const Vector& estimatedNormal,
279 double estimatedSliceThickness) const
278 { 280 {
279 bool isOpposite = false; 281 bool isOpposite = false;
280 282
281 if (points_.empty() || 283 if (points_.empty())
282 !hasSlice_ ||
283 !GeometryToolbox::IsParallelOrOpposite(isOpposite, slice.GetNormal(), geometry_.GetNormal()))
284 { 284 {
285 return false; 285 return false;
286 } 286 }
287 287 else if (hasSlice_)
288 double d = GeometryToolbox::ProjectAlongNormal(slice.GetOrigin(), geometry_.GetNormal()); 288 {
289 289 // Use the actual geometry of this specific slice
290 return (LinearAlgebra::IsNear(d, projectionAlongNormal_, 290 if (!GeometryToolbox::IsParallelOrOpposite(isOpposite, slice.GetNormal(), geometry_.GetNormal()))
291 sliceThickness_ / 2.0)); 291 {
292 return false;
293 }
294 else
295 {
296 double d = GeometryToolbox::ProjectAlongNormal(slice.GetOrigin(), geometry_.GetNormal());
297 return (LinearAlgebra::IsNear(d, projectionAlongNormal_, sliceThickness_ / 2.0));
298 }
299 }
300 else
301 {
302 // Use the estimated geometry for the global RT-STRUCT volume
303 if (!GeometryToolbox::IsParallelOrOpposite(isOpposite, slice.GetNormal(), estimatedNormal))
304 {
305 return false;
306 }
307 else
308 {
309 double d1 = GeometryToolbox::ProjectAlongNormal(slice.GetOrigin(), estimatedNormal);
310 double d2 = GeometryToolbox::ProjectAlongNormal(points_.front(), estimatedNormal);
311 return (LinearAlgebra::IsNear(d1, d2, estimatedSliceThickness / 2.0));
312 }
313 }
292 } 314 }
293 315
294 bool DicomStructureSet::Polygon::Project(double& x1, 316 bool DicomStructureSet::Polygon::Project(double& x1,
295 double& y1, 317 double& y1,
296 double& x2, 318 double& x2,
297 double& y2, 319 double& y2,
298 const CoordinateSystem3D& slice) const 320 const CoordinateSystem3D& slice,
321 const Vector& estimatedNormal,
322 double estimatedSliceThickness) const
299 { 323 {
300 if (!hasSlice_ || 324 if (!hasSlice_ ||
301 points_.size() <= 1) 325 points_.size() <= 1)
302 { 326 {
303 return false; 327 return false;
839 for (Polygons::const_iterator polygon = structure.polygons_.begin(); 863 for (Polygons::const_iterator polygon = structure.polygons_.begin();
840 polygon != structure.polygons_.end(); ++polygon) 864 polygon != structure.polygons_.end(); ++polygon)
841 { 865 {
842 const Points& points = polygon->GetPoints(); 866 const Points& points = polygon->GetPoints();
843 867
844 if (polygon->IsOnSlice(slice) && 868 if (polygon->IsOnSlice(slice, GetEstimatedNormal(), GetEstimatedSliceThickness()) &&
845 !points.empty()) 869 !points.empty())
846 { 870 {
847 chains.push_back(std::vector<ScenePoint2D>()); 871 chains.push_back(std::vector<ScenePoint2D>());
848 chains.back().reserve(points.size() + 1); 872 chains.back().reserve(points.size() + 1);
849 873
936 for (Polygons::const_iterator polygon = structure.polygons_.begin(); 960 for (Polygons::const_iterator polygon = structure.polygons_.begin();
937 polygon != structure.polygons_.end(); ++polygon) 961 polygon != structure.polygons_.end(); ++polygon)
938 { 962 {
939 double x1, y1, x2, y2; 963 double x1, y1, x2, y2;
940 964
941 if (polygon->Project(x1, y1, x2, y2, slice)) 965 if (polygon->Project(x1, y1, x2, y2, slice, GetEstimatedNormal(), GetEstimatedSliceThickness()))
942 { 966 {
943 double curZ = polygon->GetGeometryOrigin()[2]; 967 double curZ = polygon->GetGeometryOrigin()[2];
944 968
945 // x1,y1 and x2,y2 are in "slice" coordinates (the cutting plane 969 // x1,y1 and x2,y2 are in "slice" coordinates (the cutting plane
946 // geometry) 970 // geometry)