Mercurial > hg > orthanc-stl
comparison Sources/Plugin.cpp @ 28:410003c50b17
improved computation of RT-STRUCT geometry
author | Sebastien Jodogne <s.jodogne@gmail.com> |
---|---|
date | Thu, 04 Apr 2024 10:10:40 +0200 |
parents | 5fc3ea6acbd5 |
children | 62abf3c523f9 |
comparison
equal
deleted
inserted
replaced
27:bd269cfb0da7 | 28:410003c50b17 |
---|---|
24 | 24 |
25 #include "../Resources/Orthanc/Plugins/OrthancPluginCppWrapper.h" | 25 #include "../Resources/Orthanc/Plugins/OrthancPluginCppWrapper.h" |
26 | 26 |
27 #include <EmbeddedResources.h> | 27 #include <EmbeddedResources.h> |
28 | 28 |
29 #include <DicomFormat/DicomInstanceHasher.h> | |
29 #include <Compression/GzipCompressor.h> | 30 #include <Compression/GzipCompressor.h> |
30 #include <ChunkedBuffer.h> | 31 #include <ChunkedBuffer.h> |
31 #include <DicomParsing/FromDcmtkBridge.h> | 32 #include <DicomParsing/FromDcmtkBridge.h> |
32 #include <DicomParsing/ParsedDicomFile.h> | 33 #include <DicomParsing/ParsedDicomFile.h> |
33 #include <Images/ImageProcessing.h> | 34 #include <Images/ImageProcessing.h> |
50 | 51 |
51 #include <nifti1_io.h> | 52 #include <nifti1_io.h> |
52 | 53 |
53 #define ORTHANC_PLUGIN_NAME "stl" | 54 #define ORTHANC_PLUGIN_NAME "stl" |
54 | 55 |
56 static const char* const STL_SOP_CLASS_UID = "1.2.840.10008.5.1.4.1.1.104.3"; | |
57 | |
55 | 58 |
56 // Forward declaration | 59 // Forward declaration |
57 void ReadStaticAsset(std::string& target, | 60 void ReadStaticAsset(std::string& target, |
58 const std::string& path); | 61 const std::string& path); |
59 | 62 |
379 | 382 |
380 static bool MyParseDouble(double& value, | 383 static bool MyParseDouble(double& value, |
381 const std::string& s) | 384 const std::string& s) |
382 { | 385 { |
383 #if 1 | 386 #if 1 |
387 // This version is much faster, as "ParseDouble()" internally uses | |
388 // "boost::lexical_cast<double>()" | |
384 char* end = NULL; | 389 char* end = NULL; |
385 value = strtod(s.c_str(), &end); | 390 value = strtod(s.c_str(), &end); |
386 return (end == s.c_str() + s.size()); | 391 return (end == s.c_str() + s.size()); |
387 #else | 392 #else |
388 return Orthanc::SerializationToolbox::ParseDouble(value, s); | 393 return Orthanc::SerializationToolbox::ParseDouble(value, s); |
588 | 593 |
589 class StructureSet : public boost::noncopyable | 594 class StructureSet : public boost::noncopyable |
590 { | 595 { |
591 private: | 596 private: |
592 std::vector<StructurePolygon*> polygons_; | 597 std::vector<StructurePolygon*> polygons_; |
593 bool hasGeometry_; | |
594 Vector3D slicesNormal_; | |
595 double slicesSpacing_; | |
596 double minProjectionAlongNormal_; | |
597 double maxProjectionAlongNormal_; | |
598 std::string patientId_; | 598 std::string patientId_; |
599 std::string studyInstanceUid_; | 599 std::string studyInstanceUid_; |
600 std::string seriesInstanceUid_; | 600 std::string seriesInstanceUid_; |
601 std::string sopInstanceUid_; | 601 std::string sopInstanceUid_; |
602 bool hasFrameOfReferenceUid_; | 602 bool hasFrameOfReferenceUid_; |
603 std::string frameOfReferenceUid_; | 603 std::string frameOfReferenceUid_; |
604 | 604 |
605 void ComputeGeometry() | |
606 { | |
607 std::list<double> positionsList; | |
608 hasGeometry_ = false; | |
609 | |
610 for (size_t i = 0; i < polygons_.size(); i++) | |
611 { | |
612 assert(polygons_[i] != NULL); | |
613 | |
614 Vector3D n; | |
615 if (polygons_[i]->IsCoplanar(n)) | |
616 { | |
617 const Vector3D& point = polygons_[i]->GetPoint(0); | |
618 double z = Vector3D::DotProduct(point, n); | |
619 | |
620 if (!hasGeometry_) | |
621 { | |
622 hasGeometry_ = true; | |
623 slicesNormal_ = n; | |
624 minProjectionAlongNormal_ = z; | |
625 maxProjectionAlongNormal_ = z; | |
626 } | |
627 else if (!IsNear(std::abs(Vector3D::DotProduct(n, slicesNormal_)), 1)) | |
628 { | |
629 hasGeometry_ = false; | |
630 | |
631 // RT-STRUCT with non-parallel slices | |
632 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented); | |
633 } | |
634 else | |
635 { | |
636 minProjectionAlongNormal_ = std::min(minProjectionAlongNormal_, z); | |
637 maxProjectionAlongNormal_ = std::max(maxProjectionAlongNormal_, z); | |
638 } | |
639 | |
640 positionsList.push_back(Vector3D::DotProduct(n, point)); | |
641 } | |
642 } | |
643 | |
644 if (hasGeometry_) | |
645 { | |
646 std::vector<double> positions(positionsList.begin(), positionsList.end()); | |
647 assert(!positions.empty()); | |
648 | |
649 std::sort(positions.begin(), positions.end()); | |
650 RemoveDuplicateValues(positions); | |
651 assert(!positions.empty()); | |
652 | |
653 if (positions.size() == 1) | |
654 { | |
655 hasGeometry_ = false; | |
656 return; | |
657 } | |
658 | |
659 std::vector<double> offsets; | |
660 offsets.resize(positions.size() - 1); | |
661 | |
662 for (size_t i = 0; i < offsets.size(); i++) | |
663 { | |
664 offsets[i] = positions[i + 1] - positions[i]; | |
665 assert(offsets[i] > 0); | |
666 } | |
667 | |
668 std::sort(offsets.begin(), offsets.end()); | |
669 RemoveDuplicateValues(offsets); | |
670 | |
671 slicesSpacing_ = offsets[0]; | |
672 | |
673 for (size_t i = 1; i < offsets.size(); i++) | |
674 { | |
675 double d = offsets[i] / slicesSpacing_; | |
676 if (!IsNear(d, round(d))) | |
677 { | |
678 // Irregular spacing between the slices | |
679 hasGeometry_ = false; | |
680 break; | |
681 } | |
682 } | |
683 } | |
684 } | |
685 | |
686 void CheckHasGeometry() const | |
687 { | |
688 if (!hasGeometry_) | |
689 { | |
690 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls); | |
691 } | |
692 } | |
693 | |
694 public: | 605 public: |
695 explicit StructureSet(Orthanc::ParsedDicomFile& dicom) : | 606 explicit StructureSet(Orthanc::ParsedDicomFile& dicom) : |
696 hasGeometry_(false), | |
697 slicesSpacing_(0), | |
698 minProjectionAlongNormal_(0), | |
699 maxProjectionAlongNormal_(0), | |
700 hasFrameOfReferenceUid_(false) | 607 hasFrameOfReferenceUid_(false) |
701 { | 608 { |
702 DcmDataset& dataset = *dicom.GetDcmtkObject().getDataset(); | 609 DcmDataset& dataset = *dicom.GetDcmtkObject().getDataset(); |
703 patientId_ = GetStringValue(dataset, DCM_PatientID); | 610 patientId_ = GetStringValue(dataset, DCM_PatientID); |
704 studyInstanceUid_ = GetStringValue(dataset, DCM_StudyInstanceUID); | 611 studyInstanceUid_ = GetStringValue(dataset, DCM_StudyInstanceUID); |
758 polygons_[pos] = new StructurePolygon(dicom, i, j); | 665 polygons_[pos] = new StructurePolygon(dicom, i, j); |
759 } | 666 } |
760 } | 667 } |
761 | 668 |
762 assert(pos == countPolygons); | 669 assert(pos == countPolygons); |
763 | |
764 ComputeGeometry(); | |
765 } | 670 } |
766 | 671 |
767 ~StructureSet() | 672 ~StructureSet() |
768 { | 673 { |
769 for (size_t i = 0; i < polygons_.size(); i++) | 674 for (size_t i = 0; i < polygons_.size(); i++) |
793 return sopInstanceUid_; | 698 return sopInstanceUid_; |
794 } | 699 } |
795 | 700 |
796 std::string HashStudy() const | 701 std::string HashStudy() const |
797 { | 702 { |
798 std::string s; | 703 Orthanc::DicomInstanceHasher hasher(patientId_, studyInstanceUid_, seriesInstanceUid_, sopInstanceUid_); |
799 Orthanc::Toolbox::ComputeSHA1(s, patientId_ + "|" + studyInstanceUid_); | 704 return hasher.HashStudy(); |
800 return s; | |
801 } | 705 } |
802 | 706 |
803 size_t GetPolygonsCount() const | 707 size_t GetPolygonsCount() const |
804 { | 708 { |
805 return polygons_.size(); | 709 return polygons_.size(); |
816 assert(polygons_[i] != NULL); | 720 assert(polygons_[i] != NULL); |
817 return *polygons_[i]; | 721 return *polygons_[i]; |
818 } | 722 } |
819 } | 723 } |
820 | 724 |
821 bool HasGeometry() const | 725 bool HasFrameOfReferenceUid() const |
822 { | 726 { |
823 return hasGeometry_; | 727 return hasFrameOfReferenceUid_; |
728 } | |
729 | |
730 const std::string& GetFrameOfReferenceUid() const | |
731 { | |
732 if (hasFrameOfReferenceUid_) | |
733 { | |
734 return frameOfReferenceUid_; | |
735 } | |
736 else | |
737 { | |
738 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls); | |
739 } | |
740 } | |
741 }; | |
742 | |
743 | |
744 class StructureSetGeometry : public boost::noncopyable | |
745 { | |
746 private: | |
747 bool strict_; | |
748 Vector3D slicesNormal_; | |
749 double slicesSpacing_; | |
750 double minProjectionAlongNormal_; | |
751 double maxProjectionAlongNormal_; | |
752 size_t slicesCount_; | |
753 | |
754 static bool AreNormalsParallel(const Vector3D& a, | |
755 const Vector3D& b) | |
756 { | |
757 return IsNear(std::abs(Vector3D::DotProduct(a, b)), 1); | |
758 } | |
759 | |
760 bool LookupProjectionIndex(size_t& index, | |
761 double z) const | |
762 { | |
763 if (slicesCount_ == 0) | |
764 { | |
765 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | |
766 } | |
767 | |
768 if (z < minProjectionAlongNormal_ || | |
769 z > maxProjectionAlongNormal_) | |
770 { | |
771 throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange); | |
772 } | |
773 | |
774 assert(slicesSpacing_ > 0 && | |
775 minProjectionAlongNormal_ < maxProjectionAlongNormal_); | |
776 | |
777 double d = (z - minProjectionAlongNormal_) / slicesSpacing_; | |
778 | |
779 if (IsNear(d, round(d))) | |
780 { | |
781 if (d < 0.0 || | |
782 d > static_cast<double>(slicesCount_) - 1.0) | |
783 { | |
784 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | |
785 } | |
786 else | |
787 { | |
788 index = static_cast<size_t>(round(d)); | |
789 return true; | |
790 } | |
791 } | |
792 else | |
793 { | |
794 return false; | |
795 } | |
796 } | |
797 | |
798 public: | |
799 StructureSetGeometry(const StructureSet& structures, | |
800 bool strict) : | |
801 strict_(strict) | |
802 { | |
803 bool isValid = false; | |
804 | |
805 std::vector<double> projections; | |
806 projections.reserve(structures.GetPolygonsCount()); | |
807 | |
808 for (size_t i = 0; i < structures.GetPolygonsCount(); i++) | |
809 { | |
810 Vector3D normal; | |
811 if (structures.GetPolygon(i).IsCoplanar(normal)) | |
812 { | |
813 // Initialize the normal of the whole volume, if need be | |
814 if (!isValid) | |
815 { | |
816 isValid = true; | |
817 slicesNormal_ = normal; | |
818 } | |
819 | |
820 if (AreNormalsParallel(normal, slicesNormal_)) | |
821 { | |
822 // This is a valid slice (it is parallel to the normal) | |
823 const Vector3D& point = structures.GetPolygon(i).GetPoint(0); | |
824 projections.push_back(Vector3D::DotProduct(point, slicesNormal_)); | |
825 } | |
826 else | |
827 { | |
828 // RT-STRUCT with non-parallel slices | |
829 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented); | |
830 } | |
831 } | |
832 else | |
833 { | |
834 // Ignore slices that are not coplanar | |
835 } | |
836 } | |
837 | |
838 if (projections.empty()) | |
839 { | |
840 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented, | |
841 "Structure set without a valid geometry"); | |
842 } | |
843 | |
844 // Only keep unique projections | |
845 | |
846 std::sort(projections.begin(), projections.end()); | |
847 RemoveDuplicateValues(projections); | |
848 assert(!projections.empty()); | |
849 | |
850 if (projections.size() == 1) | |
851 { | |
852 // Volume with one single slice | |
853 minProjectionAlongNormal_ = projections[0]; | |
854 maxProjectionAlongNormal_ = projections[0]; | |
855 slicesSpacing_ = 1; // Arbitrary value | |
856 slicesCount_ = 1; | |
857 return; | |
858 } | |
859 | |
860 | |
861 // Compute the most probable spacing between the slices | |
862 | |
863 { | |
864 std::vector<double> spacings; | |
865 spacings.resize(projections.size() - 1); | |
866 | |
867 for (size_t i = 0; i < spacings.size(); i++) | |
868 { | |
869 spacings[i] = projections[i + 1] - projections[i]; | |
870 assert(spacings[i] > 0); | |
871 } | |
872 | |
873 std::sort(spacings.begin(), spacings.end()); | |
874 RemoveDuplicateValues(spacings); | |
875 | |
876 if (spacings.empty()) | |
877 { | |
878 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | |
879 } | |
880 | |
881 slicesSpacing_ = spacings[spacings.size() / 10]; // Take the 90% percentile of smallest spacings | |
882 assert(slicesSpacing_ > 0); | |
883 } | |
884 | |
885 | |
886 // Find the projection along the normal with the largest support | |
887 | |
888 bool first = true; | |
889 size_t bestSupport; | |
890 double bestProjection; | |
891 | |
892 std::list<size_t> candidates; | |
893 for (size_t i = 0; i < projections.size(); i++) | |
894 { | |
895 candidates.push_back(i); | |
896 } | |
897 | |
898 while (!candidates.empty()) | |
899 { | |
900 std::list<size_t> next; | |
901 | |
902 size_t countSupport = 0; | |
903 | |
904 std::list<size_t>::const_iterator it = candidates.begin(); | |
905 size_t reference = *it; | |
906 it++; | |
907 | |
908 while (it != candidates.end()) | |
909 { | |
910 double d = (projections[*it] - projections[reference]) / slicesSpacing_; | |
911 if (IsNear(d, round(d))) | |
912 { | |
913 countSupport ++; | |
914 } | |
915 else | |
916 { | |
917 next.push_back(*it); | |
918 } | |
919 | |
920 it++; | |
921 } | |
922 | |
923 if (first || | |
924 countSupport > bestSupport) | |
925 { | |
926 first = false; | |
927 bestSupport = countSupport; | |
928 bestProjection = projections[reference]; | |
929 } | |
930 | |
931 if (strict && | |
932 !next.empty()) | |
933 { | |
934 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat, | |
935 "Structure set with multiple support, which is not allowed in Strict mode"); | |
936 } | |
937 | |
938 candidates.swap(next); | |
939 } | |
940 | |
941 | |
942 // Compute the range of the projections | |
943 | |
944 minProjectionAlongNormal_ = bestProjection; | |
945 maxProjectionAlongNormal_ = bestProjection; | |
946 | |
947 for (size_t i = 0; i < projections.size(); i++) | |
948 { | |
949 double d = (projections[i] - bestProjection) / slicesSpacing_; | |
950 if (IsNear(d, round(d))) | |
951 { | |
952 minProjectionAlongNormal_ = std::min(minProjectionAlongNormal_, projections[i]); | |
953 maxProjectionAlongNormal_ = std::max(maxProjectionAlongNormal_, projections[i]); | |
954 } | |
955 } | |
956 | |
957 double d = (maxProjectionAlongNormal_ - minProjectionAlongNormal_) / slicesSpacing_; | |
958 if (IsNear(d, round(d))) | |
959 { | |
960 slicesCount_ = static_cast<size_t>(round(d)) + 1; | |
961 } | |
962 else | |
963 { | |
964 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | |
965 } | |
966 | |
967 | |
968 // Sanity check | |
969 | |
970 size_t a, b; | |
971 if (!LookupProjectionIndex(a, minProjectionAlongNormal_) || | |
972 !LookupProjectionIndex(b, maxProjectionAlongNormal_) || | |
973 a != 0 || | |
974 b + 1 != slicesCount_) | |
975 { | |
976 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | |
977 } | |
824 } | 978 } |
825 | 979 |
826 const Vector3D& GetSlicesNormal() const | 980 const Vector3D& GetSlicesNormal() const |
827 { | 981 { |
828 CheckHasGeometry(); | |
829 return slicesNormal_; | 982 return slicesNormal_; |
830 } | 983 } |
831 | 984 |
832 double GetSlicesSpacing() const | 985 double GetSlicesSpacing() const |
833 { | 986 { |
834 CheckHasGeometry(); | |
835 return slicesSpacing_; | 987 return slicesSpacing_; |
836 } | 988 } |
837 | 989 |
838 double GetMinProjectionAlongNormal() const | 990 double GetMinProjectionAlongNormal() const |
839 { | 991 { |
840 CheckHasGeometry(); | |
841 return minProjectionAlongNormal_; | 992 return minProjectionAlongNormal_; |
842 } | 993 } |
843 | 994 |
844 double GetMaxProjectionAlongNormal() const | 995 double GetMaxProjectionAlongNormal() const |
845 { | 996 { |
846 CheckHasGeometry(); | |
847 return maxProjectionAlongNormal_; | 997 return maxProjectionAlongNormal_; |
848 } | 998 } |
849 | 999 |
850 double ProjectAlongNormal(const StructurePolygon& polygon) const | 1000 bool ProjectAlongNormal(double& z, |
851 { | 1001 const StructurePolygon& polygon) const |
852 CheckHasGeometry(); | 1002 { |
853 return Vector3D::DotProduct(slicesNormal_, polygon.GetPoint(0)); | 1003 Vector3D normal; |
1004 if (polygon.IsCoplanar(normal) && | |
1005 AreNormalsParallel(normal, slicesNormal_)) | |
1006 { | |
1007 z = Vector3D::DotProduct(polygon.GetPoint(0), slicesNormal_); | |
1008 return true; | |
1009 } | |
1010 else | |
1011 { | |
1012 return false; | |
1013 } | |
854 } | 1014 } |
855 | 1015 |
856 size_t GetSlicesCount() const | 1016 size_t GetSlicesCount() const |
857 { | 1017 { |
858 CheckHasGeometry(); | 1018 return slicesCount_; |
859 double c = (maxProjectionAlongNormal_ - minProjectionAlongNormal_) / slicesSpacing_; | |
860 assert(c >= 0); | |
861 | |
862 if (!IsNear(c, round(c))) | |
863 { | |
864 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | |
865 } | |
866 else | |
867 { | |
868 return static_cast<size_t>(round(c)) + 1; | |
869 } | |
870 } | 1019 } |
871 | 1020 |
872 bool LookupSliceIndex(size_t& slice, | 1021 bool LookupSliceIndex(size_t& slice, |
873 const StructurePolygon& polygon) const | 1022 const StructurePolygon& polygon) const |
874 { | 1023 { |
875 CheckHasGeometry(); | 1024 double z; |
876 | 1025 return (ProjectAlongNormal(z, polygon) && |
877 double z = ProjectAlongNormal(polygon); | 1026 LookupProjectionIndex(slice, z)); |
878 | |
879 if (z < minProjectionAlongNormal_ || | |
880 z > maxProjectionAlongNormal_) | |
881 { | |
882 return false; | |
883 } | |
884 else | |
885 { | |
886 double c = (z - minProjectionAlongNormal_) / slicesSpacing_; | |
887 | |
888 if (IsNear(c, round(c))) | |
889 { | |
890 slice = static_cast<size_t>(round(c)); | |
891 return true; | |
892 } | |
893 else | |
894 { | |
895 return false; | |
896 } | |
897 } | |
898 } | |
899 | |
900 bool LookupReferencedSopInstanceUid(std::string& sopInstanceUid) const | |
901 { | |
902 if (HasGeometry()) | |
903 { | |
904 for (size_t i = 0; i < polygons_.size(); i++) | |
905 { | |
906 assert(polygons_[i] != NULL); | |
907 | |
908 Vector3D n; | |
909 if (polygons_[i]->IsCoplanar(n) && | |
910 Vector3D::DotProduct(n, slicesNormal_)) | |
911 { | |
912 sopInstanceUid = polygons_[i]->GetReferencedSopInstanceUid(); | |
913 return true; | |
914 } | |
915 } | |
916 } | |
917 | |
918 return false; | |
919 } | |
920 | |
921 bool HasFrameOfReferenceUid() const | |
922 { | |
923 return hasFrameOfReferenceUid_; | |
924 } | |
925 | |
926 const std::string& GetFrameOfReferenceUid() const | |
927 { | |
928 if (hasFrameOfReferenceUid_) | |
929 { | |
930 return frameOfReferenceUid_; | |
931 } | |
932 else | |
933 { | |
934 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls); | |
935 } | |
936 } | 1027 } |
937 }; | 1028 }; |
938 | 1029 |
939 | 1030 |
940 | 1031 |
1079 } | 1170 } |
1080 | 1171 |
1081 | 1172 |
1082 bool EncodeStructureSetMesh(std::string& stl, | 1173 bool EncodeStructureSetMesh(std::string& stl, |
1083 const StructureSet& structureSet, | 1174 const StructureSet& structureSet, |
1175 const StructureSetGeometry& geometry, | |
1084 const std::set<std::string>& roiNames, | 1176 const std::set<std::string>& roiNames, |
1085 unsigned int resolution, | 1177 unsigned int resolution, |
1086 bool smooth) | 1178 bool smooth) |
1087 { | 1179 { |
1088 if (!structureSet.HasGeometry()) | |
1089 { | |
1090 return false; | |
1091 } | |
1092 | |
1093 if (resolution < 1) | 1180 if (resolution < 1) |
1094 { | 1181 { |
1095 throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange); | 1182 throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange); |
1096 } | 1183 } |
1097 | 1184 |
1098 if (!IsNear(1, structureSet.GetSlicesNormal().ComputeNorm())) | 1185 if (!IsNear(1, geometry.GetSlicesNormal().ComputeNorm())) |
1099 { | 1186 { |
1100 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | 1187 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); |
1101 } | 1188 } |
1102 | 1189 |
1103 // TODO - Axes could be retrieved from the referenced CT volume | 1190 // TODO - Axes could be retrieved from the referenced CT volume |
1104 Vector3D axisX(1, 0, 0); | 1191 Vector3D axisX(1, 0, 0); |
1105 Vector3D axisY = Vector3D::CrossProduct(structureSet.GetSlicesNormal(), axisX); | 1192 Vector3D axisY = Vector3D::CrossProduct(geometry.GetSlicesNormal(), axisX); |
1106 | 1193 |
1107 if (!IsNear(1, axisX.ComputeNorm()) || | 1194 if (!IsNear(1, axisX.ComputeNorm()) || |
1108 !IsNear(1, axisY.ComputeNorm())) | 1195 !IsNear(1, axisY.ComputeNorm())) |
1109 { | 1196 { |
1110 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | 1197 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); |
1114 for (size_t i = 0; i < structureSet.GetPolygonsCount(); i++) | 1201 for (size_t i = 0; i < structureSet.GetPolygonsCount(); i++) |
1115 { | 1202 { |
1116 structureSet.GetPolygon(i).Add(extent, axisX, axisY); | 1203 structureSet.GetPolygon(i).Add(extent, axisX, axisY); |
1117 } | 1204 } |
1118 | 1205 |
1119 const int depth = structureSet.GetSlicesCount(); | 1206 const int depth = geometry.GetSlicesCount(); |
1120 | 1207 |
1121 vtkNew<vtkImageData> volume; | 1208 vtkNew<vtkImageData> volume; |
1122 volume->SetDimensions(resolution, resolution, depth); | 1209 volume->SetDimensions(resolution, resolution, depth); |
1123 volume->AllocateScalars(VTK_UNSIGNED_CHAR, 1); | 1210 volume->AllocateScalars(VTK_UNSIGNED_CHAR, 1); |
1124 | 1211 |
1133 // This polygon doesn't correspond to a ROI of interest | 1220 // This polygon doesn't correspond to a ROI of interest |
1134 continue; | 1221 continue; |
1135 } | 1222 } |
1136 | 1223 |
1137 size_t j; | 1224 size_t j; |
1138 if (!structureSet.LookupSliceIndex(j, polygon)) | 1225 if (!geometry.LookupSliceIndex(j, polygon)) |
1139 { | 1226 { |
1140 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); | 1227 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError); |
1141 } | 1228 } |
1142 | 1229 |
1143 std::vector<Orthanc::ImageProcessing::ImagePoint> points; | 1230 std::vector<Orthanc::ImageProcessing::ImagePoint> points; |
1160 } | 1247 } |
1161 | 1248 |
1162 volume->SetSpacing( | 1249 volume->SetSpacing( |
1163 extent.GetWidth() / static_cast<double>(resolution), | 1250 extent.GetWidth() / static_cast<double>(resolution), |
1164 extent.GetHeight() / static_cast<double>(resolution), | 1251 extent.GetHeight() / static_cast<double>(resolution), |
1165 (structureSet.GetMaxProjectionAlongNormal() - structureSet.GetMinProjectionAlongNormal()) / static_cast<double>(depth)); | 1252 (geometry.GetMaxProjectionAlongNormal() - geometry.GetMinProjectionAlongNormal()) / static_cast<double>(depth)); |
1166 | 1253 |
1167 // TODO | 1254 // TODO |
1168 // volume->SetOrigin() | 1255 // volume->SetOrigin() |
1169 | 1256 |
1170 return EncodeVolume(stl, volume.Get(), resolution, smooth); | 1257 return EncodeVolume(stl, volume.Get(), resolution, smooth); |
1319 { | 1406 { |
1320 static const char* const KEY_INSTANCE = "Instance"; | 1407 static const char* const KEY_INSTANCE = "Instance"; |
1321 static const char* const KEY_RESOLUTION = "Resolution"; | 1408 static const char* const KEY_RESOLUTION = "Resolution"; |
1322 static const char* const KEY_ROI_NAMES = "RoiNames"; | 1409 static const char* const KEY_ROI_NAMES = "RoiNames"; |
1323 static const char* const KEY_SMOOTH = "Smooth"; | 1410 static const char* const KEY_SMOOTH = "Smooth"; |
1411 static const char* const KEY_STRICT = "Strict"; | |
1324 | 1412 |
1325 if (request->method != OrthancPluginHttpMethod_Post) | 1413 if (request->method != OrthancPluginHttpMethod_Post) |
1326 { | 1414 { |
1327 OrthancPluginSendMethodNotAllowed(OrthancPlugins::GetGlobalContext(), output, "POST"); | 1415 OrthancPluginSendMethodNotAllowed(OrthancPlugins::GetGlobalContext(), output, "POST"); |
1328 return; | 1416 return; |
1339 Orthanc::SerializationToolbox::ReadBoolean(body, KEY_SMOOTH) : | 1427 Orthanc::SerializationToolbox::ReadBoolean(body, KEY_SMOOTH) : |
1340 true /* smooth by default */); | 1428 true /* smooth by default */); |
1341 const unsigned int resolution = (body.isMember(KEY_RESOLUTION) ? | 1429 const unsigned int resolution = (body.isMember(KEY_RESOLUTION) ? |
1342 Orthanc::SerializationToolbox::ReadUnsignedInteger(body, KEY_RESOLUTION) : | 1430 Orthanc::SerializationToolbox::ReadUnsignedInteger(body, KEY_RESOLUTION) : |
1343 256 /* default value */); | 1431 256 /* default value */); |
1432 const bool strict = (body.isMember(KEY_STRICT) ? | |
1433 Orthanc::SerializationToolbox::ReadBoolean(body, KEY_STRICT) : | |
1434 true /* strict by default */); | |
1344 | 1435 |
1345 std::set<std::string> roiNames; | 1436 std::set<std::string> roiNames; |
1346 Orthanc::SerializationToolbox::ReadSetOfStrings(roiNames, body, KEY_ROI_NAMES); | 1437 Orthanc::SerializationToolbox::ReadSetOfStrings(roiNames, body, KEY_ROI_NAMES); |
1347 | 1438 |
1348 std::unique_ptr<Orthanc::ParsedDicomFile> dicom(LoadInstance(instanceId)); | 1439 std::unique_ptr<Orthanc::ParsedDicomFile> dicom(LoadInstance(instanceId)); |
1349 | 1440 |
1350 StructureSet structureSet(*dicom); | 1441 StructureSet structureSet(*dicom); |
1351 | 1442 |
1443 StructureSetGeometry geometry(structureSet, strict); | |
1444 | |
1352 std::string stl; | 1445 std::string stl; |
1353 if (!EncodeStructureSetMesh(stl, structureSet, roiNames, resolution, smooth)) | 1446 if (!EncodeStructureSetMesh(stl, structureSet, geometry, roiNames, resolution, smooth)) |
1354 { | 1447 { |
1355 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat, "Cannot encode STL from RT-STRUCT"); | 1448 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat, "Cannot encode STL from RT-STRUCT"); |
1356 } | 1449 } |
1357 else | 1450 else |
1358 { | 1451 { |
1417 std::unique_ptr<Orthanc::ParsedDicomFile> dicom(LoadInstance(instanceId)); | 1510 std::unique_ptr<Orthanc::ParsedDicomFile> dicom(LoadInstance(instanceId)); |
1418 DcmDataset& dataset = *dicom->GetDcmtkObject().getDataset(); | 1511 DcmDataset& dataset = *dicom->GetDcmtkObject().getDataset(); |
1419 | 1512 |
1420 std::string stl; | 1513 std::string stl; |
1421 if (GetStringValue(dataset, DCM_MIMETypeOfEncapsulatedDocument) != Orthanc::MIME_STL || | 1514 if (GetStringValue(dataset, DCM_MIMETypeOfEncapsulatedDocument) != Orthanc::MIME_STL || |
1422 GetStringValue(dataset, DCM_SOPClassUID) != "1.2.840.10008.5.1.4.1.1.104.3" || | 1515 GetStringValue(dataset, DCM_SOPClassUID) != STL_SOP_CLASS_UID || |
1423 !dicom->GetTagValue(stl, Orthanc::DICOM_TAG_ENCAPSULATED_DOCUMENT)) | 1516 !dicom->GetTagValue(stl, Orthanc::DICOM_TAG_ENCAPSULATED_DOCUMENT)) |
1424 { | 1517 { |
1425 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadRequest, "DICOM instance not encapsulating a STL model: " + instanceId); | 1518 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadRequest, "DICOM instance not encapsulating a STL model: " + instanceId); |
1426 } | 1519 } |
1427 else | 1520 else |