orthanc-stone: OrthancStone/Sources/Toolbox/DicomStructureSet.cpp comparison

comparison OrthancStone/Sources/Toolbox/DicomStructureSet.cpp @ 2161:e65fe2e50fde dicom-sr tip

integration mainline->dicom-sr

author	Sebastien Jodogne <s.jodogne@gmail.com>
date	Fri, 27 Sep 2024 22:34:17 +0200
parents	917e40af6b45
children

comparison

equal deleted inserted replaced

-:f68f9a8d0d63
+:e65fe2e50fde
 #endif
 namespace OrthancStone
 {
+static const Orthanc::DicomTag DICOM_TAG_CONTOUR_DATA(0x3006, 0x0050);
 static const Orthanc::DicomTag DICOM_TAG_CONTOUR_GEOMETRIC_TYPE(0x3006, 0x0042);
 static const Orthanc::DicomTag DICOM_TAG_CONTOUR_IMAGE_SEQUENCE(0x3006, 0x0016);
 static const Orthanc::DicomTag DICOM_TAG_CONTOUR_SEQUENCE(0x3006, 0x0040);
-static const Orthanc::DicomTag DICOM_TAG_CONTOUR_DATA(0x3006, 0x0050);
 static const Orthanc::DicomTag DICOM_TAG_NUMBER_OF_CONTOUR_POINTS(0x3006, 0x0046);
+static const Orthanc::DicomTag DICOM_TAG_REFERENCED_ROI_NUMBER(0x3006, 0x0084);
 static const Orthanc::DicomTag DICOM_TAG_REFERENCED_SOP_INSTANCE_UID(0x0008, 0x1155);
 static const Orthanc::DicomTag DICOM_TAG_ROI_CONTOUR_SEQUENCE(0x3006, 0x0039);
 static const Orthanc::DicomTag DICOM_TAG_ROI_DISPLAY_COLOR(0x3006, 0x002a);
 static const Orthanc::DicomTag DICOM_TAG_ROI_NAME(0x3006, 0x0026);
+static const Orthanc::DicomTag DICOM_TAG_ROI_NUMBER(0x3006, 0x0022);
 static const Orthanc::DicomTag DICOM_TAG_RT_ROI_INTERPRETED_TYPE(0x3006, 0x00a4);
 static const Orthanc::DicomTag DICOM_TAG_RT_ROI_OBSERVATIONS_SEQUENCE(0x3006, 0x0080);
 static const Orthanc::DicomTag DICOM_TAG_STRUCTURE_SET_ROI_SEQUENCE(0x3006, 0x0020);
 }
 }
 }
-static void AddSegmentIfIntersection(Extent2D& extent,
+void DicomStructureSet::Polygon::Project(std::list<Extent2D>& target,
-const CoordinateSystem3D& cuttingPlane,
-const Vector& p1,
-const Vector& p2,
-double originDistance)
-{
-// Does this segment intersects the cutting plane?
-double d1 = cuttingPlane.ProjectAlongNormal(p1);
-double d2 = cuttingPlane.ProjectAlongNormal(p2);
-if ((d1 < originDistance && d2 > originDistance) ||
-(d1 > originDistance && d2 < originDistance))
-{
-// This is an intersection: Add the segment
-double x, y;
-cuttingPlane.ProjectPoint2(x, y, p1);
-extent.AddPoint(x, y);
-cuttingPlane.ProjectPoint2(x, y, p2);
-extent.AddPoint(x, y);
-}
-}
-bool DicomStructureSet::Polygon::Project(double& x1,
-double& y1,
-double& x2,
-double& y2,
 const CoordinateSystem3D& cuttingPlane,
 const Vector& estimatedNormal,
 double estimatedSliceThickness) const
 {
-if (points_.size() <= 1)
+CoordinateSystem3D geometry;
-{
-return false;
-}
-Vector normal = estimatedNormal;
 double thickness = estimatedSliceThickness;
+/**
+* 1. Estimate the 3D plane associated with this polygon.
+**/
 if (hasSlice_)
 {
-normal = geometry_.GetNormal();
+// The exact geometry is known for this slice
+geometry = geometry_;
 thickness = sliceThickness_;
 }
+else if (points_.size() < 2)
+{
+return;
+}
+else
+{
+// Estimate the geometry
+Vector origin = points_[0];
+Vector axisX;
+bool found = false;
+for (size_t i = 1; i < points_.size(); i++)
+{
+axisX = points_[1] - origin;
+if (boost::numeric::ublas::norm_2(axisX) > 10.0 * std::numeric_limits<double>::epsilon())
+{
+found = true;
+break;
+}
+}
+if (!found)
+{
+return;  // The polygon is too small to extract a reliable geometry out of it
+}
+LinearAlgebra::NormalizeVector(axisX);
+Vector axisY;
+LinearAlgebra::CrossProduct(axisY, axisX, estimatedNormal);
+geometry = CoordinateSystem3D(origin, axisX, axisY);
+}
+/**
+* 2. Project the 3D cutting plane as a 2D line onto the polygon plane.
+**/
+double cuttingX1, cuttingY1, cuttingX2, cuttingY2;
+geometry.ProjectPoint(cuttingX1, cuttingY1, cuttingPlane.GetOrigin());
 bool isOpposite;
-if (!GeometryToolbox::IsParallelOrOpposite(isOpposite, normal, cuttingPlane.GetAxisX()) &&
+if (GeometryToolbox::IsParallelOrOpposite(isOpposite, geometry.GetNormal(), cuttingPlane.GetAxisX()))
-!GeometryToolbox::IsParallelOrOpposite(isOpposite, normal, cuttingPlane.GetAxisY()))
+{
-{
+geometry.ProjectPoint(cuttingX2, cuttingY2, cuttingPlane.GetOrigin() + cuttingPlane.GetAxisY());
-return false;
+}
-}
+else if (GeometryToolbox::IsParallelOrOpposite(isOpposite, geometry.GetNormal(), cuttingPlane.GetAxisY()))
+{
-const double d = cuttingPlane.ProjectAlongNormal(cuttingPlane.GetOrigin());
+geometry.ProjectPoint(cuttingX2, cuttingY2, cuttingPlane.GetOrigin() + cuttingPlane.GetAxisX());
-Extent2D extent;
-AddSegmentIfIntersection(extent, cuttingPlane, points_[points_.size() - 1], points_[0], d);
-for (size_t i = 1; i < points_.size(); i++)
-{
-AddSegmentIfIntersection(extent, cuttingPlane, points_[i - 1], points_[i], d);
-}
-if (extent.GetWidth() > 0 ||
-extent.GetHeight() > 0)
-{
-// Let's convert them to 3D world geometry to add the slice thickness
-Vector p1 = (cuttingPlane.MapSliceToWorldCoordinates(extent.GetX1(), extent.GetY1()) +
-thickness / 2.0 * normal);
-Vector p2 = (cuttingPlane.MapSliceToWorldCoordinates(extent.GetX2(), extent.GetY2()) -
-thickness / 2.0 * normal);
-// Then back to the cutting plane geometry
-cuttingPlane.ProjectPoint2(x1, y1, p1);
-cuttingPlane.ProjectPoint2(x2, y2, p2);
-return true;
 }
 else
 {
-return false;
+return;
+}
+/**
+* 3. Compute the intersection of the 2D cutting line with the polygon.
+**/
+// Initialize the projection of a point onto a line:
+// https://stackoverflow.com/a/64330724
+const double abx = cuttingX2 - cuttingX1;
+const double aby = cuttingY2 - cuttingY1;
+const double denominator = abx * abx + aby * aby;
+if (LinearAlgebra::IsCloseToZero(denominator))
+{
+return;  // Should never happen
+}
+std::vector<double> intersections;
+intersections.reserve(points_.size());
+for (size_t i = 0; i < points_.size(); i++)
+{
+double segmentX1, segmentY1, segmentX2, segmentY2;
+geometry.ProjectPoint(segmentX1, segmentY1, points_[i]);
+geometry.ProjectPoint(segmentX2, segmentY2, points_[(i + 1) % points_.size()]);
+double x, y;
+if (GeometryToolbox::IntersectLineAndSegment(x, y, cuttingX1, cuttingY1, cuttingX2, cuttingY2,
+segmentX1, segmentY1, segmentX2, segmentY2))
+{
+// For each polygon segment that intersect the cutting line,
+// register its offset over the cutting line
+const double acx = x - cuttingX1;
+const double acy = y - cuttingY1;
+intersections.push_back((abx * acx + aby * acy) / denominator);
+}
+}
+/**
+* 4. Sort the intersection offsets, then generates one 2D rectangle on the
+* cutting plane from each pair of successive intersections.
+**/
+std::sort(intersections.begin(), intersections.end());
+if (intersections.size() % 2 == 1)
+{
+return;  // Should never happen
+}
+for (size_t i = 0; i < intersections.size(); i += 2)
+{
+Vector p1, p2;
+{
+// Let's convert them to 3D world geometry to add the slice thickness
+const double x1 = cuttingX1 + intersections[i] * abx;
+const double y1 = cuttingY1 + intersections[i] * aby;
+const double x2 = cuttingX1 + intersections[i + 1] * abx;
+const double y2 = cuttingY1 + intersections[i + 1] * aby;
+p1 = (geometry.MapSliceToWorldCoordinates(x1, y1) + thickness / 2.0 * geometry.GetNormal());
+p2 = (geometry.MapSliceToWorldCoordinates(x2, y2) - thickness / 2.0 * geometry.GetNormal());
+}
+{
+// Then back to the cutting plane geometry
+double x1, y1, x2, y2;
+cuttingPlane.ProjectPoint2(x1, y1, p1);
+cuttingPlane.ProjectPoint2(x2, y2, p2);
+target.push_back(Extent2D(x1, y1, x2, y2));
+}
 }
 }
 const DicomStructureSet::Structure& DicomStructureSet::GetStructure(size_t index) const
 {
 #if STONE_TIME_BLOCKING_OPS
 boost::posix_time::ptime timerStart = boost::posix_time::microsec_clock::universal_time();
 #endif
+std::map<int, size_t> roiNumbersIndex;
 DicomDatasetReader reader(tags);
-size_t count, tmp;
-if (!tags.GetSequenceSize(count, Orthanc::DicomPath(DICOM_TAG_RT_ROI_OBSERVATIONS_SEQUENCE)) ||
+/**
-!tags.GetSequenceSize(tmp, Orthanc::DicomPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE)) ||
+* 1. Read all the available ROIs.
-tmp != count ||
+**/
-!tags.GetSequenceSize(tmp, Orthanc::DicomPath(DICOM_TAG_STRUCTURE_SET_ROI_SEQUENCE)) ||
-tmp != count)
+{
-{
+size_t count;
-throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+if (!tags.GetSequenceSize(count, Orthanc::DicomPath(DICOM_TAG_STRUCTURE_SET_ROI_SEQUENCE)))
-}
+{
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
-structures_.resize(count);
+}
-structureNamesIndex_.clear();
+structures_.resize(count);
-for (size_t i = 0; i < count; i++)
+structureNamesIndex_.clear();
-{
-structures_[i].interpretation_ = reader.GetStringValue
+for (size_t i = 0; i < count; i++)
-(Orthanc::DicomPath(DICOM_TAG_RT_ROI_OBSERVATIONS_SEQUENCE, i,
+{
-DICOM_TAG_RT_ROI_INTERPRETED_TYPE),
+int roiNumber;
-"No interpretation");
+if (!reader.GetIntegerValue
+(roiNumber, Orthanc::DicomPath(DICOM_TAG_STRUCTURE_SET_ROI_SEQUENCE, i, DICOM_TAG_ROI_NUMBER)))
-structures_[i].name_ = reader.GetStringValue
+{
-(Orthanc::DicomPath(DICOM_TAG_STRUCTURE_SET_ROI_SEQUENCE, i,
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
-DICOM_TAG_ROI_NAME),
+}
-"No name");
+if (roiNumbersIndex.find(roiNumber) != roiNumbersIndex.end())
-if (structureNamesIndex_.find(structures_[i].name_) == structureNamesIndex_.end())
+{
-{
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat,
-structureNamesIndex_[structures_[i].name_] = i;
+"Twice the same ROI number: " + boost::lexical_cast<std::string>(roiNumber));
 }
-else
-{
+roiNumbersIndex[roiNumber] = i;
-throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat,
-"RT-STRUCT with twice the same name for a structure: " + structures_[i].name_);
+structures_[i].name_ = reader.GetStringValue
-}
+(Orthanc::DicomPath(DICOM_TAG_STRUCTURE_SET_ROI_SEQUENCE, i, DICOM_TAG_ROI_NAME), "No name");
+structures_[i].interpretation_ = "No interpretation";
-Vector color;
-if (FastParseVector(color, tags, Orthanc::DicomPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+if (structureNamesIndex_.find(structures_[i].name_) == structureNamesIndex_.end())
-DICOM_TAG_ROI_DISPLAY_COLOR)) &&
+{
-color.size() == 3)
+structureNamesIndex_[structures_[i].name_] = i;
-{
+}
-structures_[i].red_ = ConvertColor(color[0]);
+else
-structures_[i].green_ = ConvertColor(color[1]);
+{
-structures_[i].blue_ = ConvertColor(color[2]);
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat,
-}
+"RT-STRUCT with twice the same name for a structure: " + structures_[i].name_);
-else
+}
-{
+}
-structures_[i].red_ = 255;
+}
-structures_[i].green_ = 0;
-structures_[i].blue_ = 0;
-}
+/**
+* 2. Read the interpretation of the ROIs (if available).
-size_t countSlices;
+**/
-if (!tags.GetSequenceSize(countSlices, Orthanc::DicomPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
-DICOM_TAG_CONTOUR_SEQUENCE)))
+{
-{
+size_t count;
-countSlices = 0;
+if (!tags.GetSequenceSize(count, Orthanc::DicomPath(DICOM_TAG_RT_ROI_OBSERVATIONS_SEQUENCE)))
-}
+{
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
-LOG(INFO) << "New RT structure: \"" << structures_[i].name_
+}
-<< "\" with interpretation \"" << structures_[i].interpretation_
-<< "\" containing " << countSlices << " slices (color: "
+for (size_t i = 0; i < count; i++)
-<< static_cast<int>(structures_[i].red_) << ","
+{
-<< static_cast<int>(structures_[i].green_) << ","
+std::string interpretation;
-<< static_cast<int>(structures_[i].blue_) << ")";
+if (reader.GetDataset().GetStringValue(interpretation,
+Orthanc::DicomPath(DICOM_TAG_RT_ROI_OBSERVATIONS_SEQUENCE, i,
-/**
+DICOM_TAG_RT_ROI_INTERPRETED_TYPE)))
-* These temporary variables avoid allocating many vectors in
+{
-* the loop below (indeed, "Orthanc::DicomPath" handles a
+int roiNumber;
-* "std::vector<PrefixItem>")
+if (!reader.GetIntegerValue(roiNumber,
-**/
+Orthanc::DicomPath(DICOM_TAG_RT_ROI_OBSERVATIONS_SEQUENCE, i,
-Orthanc::DicomPath countPointsPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+DICOM_TAG_REFERENCED_ROI_NUMBER)))
-DICOM_TAG_CONTOUR_SEQUENCE, 0,
+{
-DICOM_TAG_NUMBER_OF_CONTOUR_POINTS);
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+}
-Orthanc::DicomPath geometricTypePath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+std::map<int, size_t>::const_iterator found = roiNumbersIndex.find(roiNumber);
+if (found == roiNumbersIndex.end())
+{
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+}
+structures_[found->second].interpretation_ = interpretation;
+}
+}
+}
+/**
+* 3. Read the contours.
+**/
+{
+size_t count;
+if (!tags.GetSequenceSize(count, Orthanc::DicomPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE)))
+{
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+}
+for (size_t i = 0; i < count; i++)
+{
+int roiNumber;
+if (!reader.GetIntegerValue(roiNumber,
+Orthanc::DicomPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+DICOM_TAG_REFERENCED_ROI_NUMBER)))
+{
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+}
+std::map<int, size_t>::const_iterator found = roiNumbersIndex.find(roiNumber);
+if (found == roiNumbersIndex.end())
+{
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+}
+Structure& target = structures_[found->second];
+Vector color;
+if (FastParseVector(color, tags, Orthanc::DicomPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+DICOM_TAG_ROI_DISPLAY_COLOR)) &&
+color.size() == 3)
+{
+target.red_ = ConvertColor(color[0]);
+target.green_ = ConvertColor(color[1]);
+target.blue_ = ConvertColor(color[2]);
+}
+else
+{
+target.red_ = 255;
+target.green_ = 0;
+target.blue_ = 0;
+}
+size_t countSlices;
+if (!tags.GetSequenceSize(countSlices, Orthanc::DicomPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+DICOM_TAG_CONTOUR_SEQUENCE)))
+{
+countSlices = 0;
+}
+LOG(INFO) << "New RT structure: \"" << target.name_
+<< "\" with interpretation \"" << target.interpretation_
+<< "\" containing " << countSlices << " slices (color: "
+<< static_cast<int>(target.red_) << ","
+<< static_cast<int>(target.green_) << ","
+<< static_cast<int>(target.blue_) << ")";
+/**
+* These temporary variables avoid allocating many vectors in
+* the loop below (indeed, "Orthanc::DicomPath" handles a
+* "std::vector<PrefixItem>")
+**/
+Orthanc::DicomPath countPointsPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
 DICOM_TAG_CONTOUR_SEQUENCE, 0,
-DICOM_TAG_CONTOUR_GEOMETRIC_TYPE);
+DICOM_TAG_NUMBER_OF_CONTOUR_POINTS);
-Orthanc::DicomPath imageSequencePath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+Orthanc::DicomPath geometricTypePath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+DICOM_TAG_CONTOUR_SEQUENCE, 0,
+DICOM_TAG_CONTOUR_GEOMETRIC_TYPE);
+Orthanc::DicomPath imageSequencePath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+DICOM_TAG_CONTOUR_SEQUENCE, 0,
+DICOM_TAG_CONTOUR_IMAGE_SEQUENCE);
+// (3006,0039)[i] / (0x3006, 0x0040)[0] / (0x3006, 0x0016)[0] / (0x0008, 0x1155)
+Orthanc::DicomPath referencedInstancePath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+DICOM_TAG_CONTOUR_SEQUENCE, 0,
+DICOM_TAG_CONTOUR_IMAGE_SEQUENCE, 0,
+DICOM_TAG_REFERENCED_SOP_INSTANCE_UID);
+Orthanc::DicomPath contourDataPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
 DICOM_TAG_CONTOUR_SEQUENCE, 0,
-DICOM_TAG_CONTOUR_IMAGE_SEQUENCE);
+DICOM_TAG_CONTOUR_DATA);
-// (3006,0039)[i] / (0x3006, 0x0040)[0] / (0x3006, 0x0016)[0] / (0x0008, 0x1155)
+for (size_t j = 0; j < countSlices; j++)
-Orthanc::DicomPath referencedInstancePath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+{
-DICOM_TAG_CONTOUR_SEQUENCE, 0,
+unsigned int countPoints;
-DICOM_TAG_CONTOUR_IMAGE_SEQUENCE, 0,
-DICOM_TAG_REFERENCED_SOP_INSTANCE_UID);
+countPointsPath.SetPrefixIndex(1, j);
+if (!reader.GetUnsignedIntegerValue(countPoints, countPointsPath))
-Orthanc::DicomPath contourDataPath(DICOM_TAG_ROI_CONTOUR_SEQUENCE, i,
+{
-DICOM_TAG_CONTOUR_SEQUENCE, 0,
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
-DICOM_TAG_CONTOUR_DATA);
+}
-for (size_t j = 0; j < countSlices; j++)
+//LOG(INFO) << "Parsing slice containing " << countPoints << " vertices";
-{
-unsigned int countPoints;
+geometricTypePath.SetPrefixIndex(1, j);
+std::string type = reader.GetMandatoryStringValue(geometricTypePath);
-countPointsPath.SetPrefixIndex(1, j);
+if (type != "CLOSED_PLANAR")
-if (!reader.GetUnsignedIntegerValue(countPoints, countPointsPath))
+{
-{
+LOG(WARNING) << "Ignoring contour with geometry type: " << type;
-throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+continue;
 }
-//LOG(INFO) << "Parsing slice containing " << countPoints << " vertices";
+size_t size;
-geometricTypePath.SetPrefixIndex(1, j);
+imageSequencePath.SetPrefixIndex(1, j);
-std::string type = reader.GetMandatoryStringValue(geometricTypePath);
+if (!tags.GetSequenceSize(size, imageSequencePath) || size != 1)
-if (type != "CLOSED_PLANAR")
+{
-{
+LOG(ERROR) << "The ContourImageSequence sequence (tag 3006,0016) must be present and contain one entry.";
-LOG(WARNING) << "Ignoring contour with geometry type: " << type;
+throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
-continue;
+}
-}
+referencedInstancePath.SetPrefixIndex(1, j);
-size_t size;
+std::string sopInstanceUid = reader.GetMandatoryStringValue(referencedInstancePath);
-imageSequencePath.SetPrefixIndex(1, j);
+contourDataPath.SetPrefixIndex(1, j);
-if (!tags.GetSequenceSize(size, imageSequencePath) || size != 1)
+std::string slicesData = reader.GetMandatoryStringValue(contourDataPath);
-{
-LOG(ERROR) << "The ContourImageSequence sequence (tag 3006,0016) must be present and contain one entry.";
+Vector points;
-throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
-}
+if (!GenericToolbox::FastParseVector(points, slicesData) ||
+points.size() != 3 * countPoints)
-referencedInstancePath.SetPrefixIndex(1, j);
+{
-std::string sopInstanceUid = reader.GetMandatoryStringValue(referencedInstancePath);
+throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+}
-contourDataPath.SetPrefixIndex(1, j);
-std::string slicesData = reader.GetMandatoryStringValue(contourDataPath);
+// seen in real world
+if(Orthanc::Toolbox::StripSpaces(sopInstanceUid) == "")
-Vector points;
+{
+LOG(ERROR) << "WARNING. The following Dicom tag (Referenced SOP Instance UID) contains an empty value : // (3006,0039)[" << i << "] / (0x3006, 0x0040)[0] / (0x3006, 0x0016)[0] / (0x0008, 0x1155)";
-if (!GenericToolbox::FastParseVector(points, slicesData) ||
+}
-points.size() != 3 * countPoints)
-{
+Polygon polygon(sopInstanceUid);
-throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
+polygon.Reserve(countPoints);
-}
+for (size_t k = 0; k < countPoints; k++)
-// seen in real world
+{
-if(Orthanc::Toolbox::StripSpaces(sopInstanceUid) == "")
+Vector v(3);
-{
+v[0] = points[3 * k];
-LOG(ERROR) << "WARNING. The following Dicom tag (Referenced SOP Instance UID) contains an empty value : // (3006,0039)[" << i << "] / (0x3006, 0x0040)[0] / (0x3006, 0x0016)[0] / (0x0008, 0x1155)";
+v[1] = points[3 * k + 1];
-}
+v[2] = points[3 * k + 2];
+polygon.AddPoint(v);
-Polygon polygon(sopInstanceUid);
+}
-polygon.Reserve(countPoints);
+target.polygons_.push_back(polygon);
-for (size_t k = 0; k < countPoints; k++)
+}
-{
-Vector v(3);
-v[0] = points[3 * k];
-v[1] = points[3 * k + 1];
-v[2] = points[3 * k + 2];
-polygon.AddPoint(v);
-}
-structures_[i].polygons_.push_back(polygon);
 }
 }
 EstimateGeometry();
 std::vector<BoostPolygon> projected;
 for (Polygons::const_iterator polygon = structure.polygons_.begin();
 polygon != structure.polygons_.end(); ++polygon)
 {
-double x1, y1, x2, y2;
+std::list<Extent2D> rectangles;
+polygon->Project(rectangles, cutting, GetEstimatedNormal(), GetEstimatedSliceThickness());
-if (polygon->Project(x1, y1, x2, y2, cutting, GetEstimatedNormal(), GetEstimatedSliceThickness()))
-{
+for (std::list<Extent2D>::const_iterator it = rectangles.begin(); it != rectangles.end(); ++it)
-projected.push_back(CreateRectangle(x1, y1, x2, y2));
+{
+projected.push_back(CreateRectangle(it->GetX1(), it->GetY1(), it->GetX2(), it->GetY2()));
 }
 }
 BoostMultiPolygon merged;
 Union(merged, projected);
 std::list<Extent2D> rectangles;
 for (Polygons::const_iterator polygon = structure.polygons_.begin();
 polygon != structure.polygons_.end(); ++polygon)
 {
-double x1, y1, x2, y2;
+polygon->Project(rectangles, cutting, GetEstimatedNormal(), GetEstimatedSliceThickness());
-if (polygon->Project(x1, y1, x2, y2, cutting, GetEstimatedNormal(), GetEstimatedSliceThickness()))
-{
-rectangles.push_back(Extent2D(x1, y1, x2, y2));
-}
 }
 typedef std::list< std::vector<ScenePoint2D> >  Contours;
 Contours contours;

Mercurial > hg > orthanc-stone

comparison OrthancStone/Sources/Toolbox/DicomStructureSet.cpp @ 2161:e65fe2e50fde dicom-sr tip