diff --git a/Modules/SegmentationUI/Qmitk/QmitkSlicesInterpolator.cpp b/Modules/SegmentationUI/Qmitk/QmitkSlicesInterpolator.cpp index 86f1b19d0c..09f0675a1c 100644 --- a/Modules/SegmentationUI/Qmitk/QmitkSlicesInterpolator.cpp +++ b/Modules/SegmentationUI/Qmitk/QmitkSlicesInterpolator.cpp @@ -1,1999 +1,2007 @@ /*============================================================================ The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center (DKFZ) All rights reserved. Use of this source code is governed by a 3-clause BSD license that can be found in the LICENSE file. ============================================================================*/ #include "QmitkSlicesInterpolator.h" #include "QmitkRenderWindow.h" #include "QmitkRenderWindowWidget.h" #include "mitkApplyDiffImageOperation.h" #include "mitkColorProperty.h" #include "mitkCoreObjectFactory.h" #include "mitkDiffImageApplier.h" #include "mitkInteractionConst.h" #include "mitkLevelWindowProperty.h" #include "mitkOperationEvent.h" #include "mitkProgressBar.h" #include "mitkProperties.h" #include "mitkRenderingManager.h" #include "mitkSegTool2D.h" #include "mitkSliceNavigationController.h" #include "mitkSurfaceToImageFilter.h" #include "mitkToolManager.h" #include "mitkUndoController.h" #include #include #include #include #include #include #include #include #include #include #include #include // Includes for the merge operation #include "mitkImageToContourFilter.h" #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include namespace { template itk::SmartPointer GetData(const mitk::DataNode* dataNode) { return nullptr != dataNode ? dynamic_cast(dataNode->GetData()) : nullptr; } } float SURFACE_COLOR_RGB[3] = {0.49f, 1.0f, 0.16f}; const std::map QmitkSlicesInterpolator::createActionToSlicer(const QList& windows) { std::map actionToSliceDimension; for (auto* window : windows) { std::string windowName; auto renderWindowWidget = dynamic_cast(window->parentWidget()); if (renderWindowWidget) { windowName = renderWindowWidget->GetCornerAnnotationText(); } else { windowName = window->GetRenderer()->GetName(); } auto slicer = window->GetSliceNavigationController(); actionToSliceDimension[new QAction(QString::fromStdString(windowName), nullptr)] = slicer; } return actionToSliceDimension; } // Check whether the given contours are coplanar bool AreContoursCoplanar(mitk::SurfaceInterpolationController::ContourPositionInformation leftHandSide, mitk::SurfaceInterpolationController::ContourPositionInformation rightHandSide) { // Here we check two things: // 1. Whether the normals of both contours are at least parallel // 2. Whether both contours lie in the same plane // Check for coplanarity: // a. Span a vector between two points one from each contour // b. Calculate dot product for the vector and one of the normals // c. If the dot is zero the two vectors are orthogonal and the contours are coplanar double vec[3]; vec[0] = leftHandSide.ContourPoint[0] - rightHandSide.ContourPoint[0]; vec[1] = leftHandSide.ContourPoint[1] - rightHandSide.ContourPoint[1]; vec[2] = leftHandSide.ContourPoint[2] - rightHandSide.ContourPoint[2]; double n[3]; n[0] = rightHandSide.ContourNormal[0]; n[1] = rightHandSide.ContourNormal[1]; n[2] = rightHandSide.ContourNormal[2]; double dot = vtkMath::Dot(n, vec); double n2[3]; n2[0] = leftHandSide.ContourNormal[0]; n2[1] = leftHandSide.ContourNormal[1]; n2[2] = leftHandSide.ContourNormal[2]; // The normals of both contours have to be parallel but not of the same orientation double lengthLHS = leftHandSide.ContourNormal.GetNorm(); double lengthRHS = rightHandSide.ContourNormal.GetNorm(); double dot2 = vtkMath::Dot(n, n2); bool contoursParallel = mitk::Equal(fabs(lengthLHS * lengthRHS), fabs(dot2), 0.001); if (mitk::Equal(dot, 0.0, 0.001) && contoursParallel) return true; else return false; } mitk::Image::Pointer ExtractSliceFromImage(mitk::Image* image, const mitk::PlaneGeometry * contourPlane, unsigned int timeStep) { vtkSmartPointer reslice = vtkSmartPointer::New(); // set to false to extract a slice reslice->SetOverwriteMode(false); reslice->Modified(); mitk::ExtractSliceFilter::Pointer extractor = mitk::ExtractSliceFilter::New(reslice); extractor->SetInput(image); extractor->SetTimeStep(timeStep); extractor->SetWorldGeometry(contourPlane); extractor->SetVtkOutputRequest(false); extractor->SetResliceTransformByGeometry(image->GetTimeGeometry()->GetGeometryForTimeStep(timeStep)); extractor->Update(); mitk::Image::Pointer slice = extractor->GetOutput(); return slice; } template std::vector GetPixelValuesPresentInImage(mitk::LabelSetImage* labelSetImage) { std::vector pixelsPresent; mitk::ImagePixelReadAccessor readAccessor(labelSetImage); std::size_t numberOfPixels = 1; for (size_t dim = 0; dim < VImageDimension; ++dim) numberOfPixels *= static_cast(readAccessor.GetDimension(dim)); auto src = readAccessor.GetData(); for (std::size_t i = 0; i < numberOfPixels; ++i) { mitk::Label::PixelType pixelVal = *(src + i); if ( (std::find(pixelsPresent.begin(), pixelsPresent.end(), pixelVal) == pixelsPresent.end()) && (pixelVal != mitk::LabelSetImage::UnlabeledValue) ) pixelsPresent.push_back(pixelVal); } return pixelsPresent; } template ModifyLabelActionTrigerred ModifyLabelProcessing(mitk::LabelSetImage* labelSetImage, mitk::SurfaceInterpolationController::Pointer surfaceInterpolator, unsigned int timePoint) { auto currentLayerID = labelSetImage->GetActiveLayer(); auto numTimeSteps = labelSetImage->GetTimeSteps(); ModifyLabelActionTrigerred actionTriggered = ModifyLabelActionTrigerred::Null; - mitk::SurfaceInterpolationController::ContourPositionInformationList ¤tContourList = - surfaceInterpolator->GetContours(timePoint, currentLayerID); + auto* currentContourList = surfaceInterpolator->GetContours(timePoint, currentLayerID); + + if (nullptr == currentContourList) + { + surfaceInterpolator->OnAddLayer(); + currentContourList = surfaceInterpolator->GetContours(timePoint, currentLayerID); + } mitk::LabelSetImage::Pointer labelSetImage2 = labelSetImage->Clone(); mitk::ImagePixelReadAccessor readAccessor(labelSetImage2.GetPointer()); - for (auto& contour : currentContourList) + for (auto& contour : *currentContourList) { mitk::Label::PixelType contourPixelValue; itk::Index<3> itkIndex; labelSetImage2->GetGeometry()->WorldToIndex(contour.ContourPoint, itkIndex); if (VImageDimension == 4) { itk::Index time3DIndex; for (size_t i = 0; i < itkIndex.size(); ++i) time3DIndex[i] = itkIndex[i]; time3DIndex[3] = timePoint; contourPixelValue = readAccessor.GetPixelByIndexSafe(time3DIndex); } else if (VImageDimension == 3) { itk::Index geomIndex; for (size_t i = 0; i < itkIndex.size(); ++i) geomIndex[i] = itkIndex[i]; contourPixelValue = readAccessor.GetPixelByIndexSafe(geomIndex); } if (contour.LabelValue != contourPixelValue) { if (contourPixelValue == 0) // Erase label { for (size_t t = 0; t < numTimeSteps; ++t) surfaceInterpolator->RemoveContours(contour.LabelValue, t, currentLayerID); actionTriggered = ModifyLabelActionTrigerred::Erase; } else { contour.LabelValue = contourPixelValue; actionTriggered = ModifyLabelActionTrigerred::Merge; } } } return actionTriggered; } QmitkSlicesInterpolator::QmitkSlicesInterpolator(QWidget *parent, const char * /*name*/) : QWidget(parent), m_Interpolator(mitk::SegmentationInterpolationController::New()), m_SurfaceInterpolator(mitk::SurfaceInterpolationController::GetInstance()), m_ToolManager(nullptr), m_Initialized(false), m_LastSNC(nullptr), m_LastSliceIndex(0), m_2DInterpolationEnabled(false), m_3DInterpolationEnabled(false), m_PreviousActiveLabelValue(0), m_CurrentActiveLabelValue(0), m_PreviousLayerIndex(0), m_CurrentLayerIndex(0), m_FirstRun(true) { m_GroupBoxEnableExclusiveInterpolationMode = new QGroupBox("Interpolation", this); QVBoxLayout *vboxLayout = new QVBoxLayout(m_GroupBoxEnableExclusiveInterpolationMode); m_EdgeDetector = mitk::FeatureBasedEdgeDetectionFilter::New(); m_PointScorer = mitk::PointCloudScoringFilter::New(); m_CmbInterpolation = new QComboBox(m_GroupBoxEnableExclusiveInterpolationMode); m_CmbInterpolation->addItem("Disabled"); m_CmbInterpolation->addItem("2-Dimensional"); m_CmbInterpolation->addItem("3-Dimensional"); vboxLayout->addWidget(m_CmbInterpolation); m_BtnApply2D = new QPushButton("Confirm for single slice", m_GroupBoxEnableExclusiveInterpolationMode); vboxLayout->addWidget(m_BtnApply2D); m_BtnApplyForAllSlices2D = new QPushButton("Confirm for all slices", m_GroupBoxEnableExclusiveInterpolationMode); vboxLayout->addWidget(m_BtnApplyForAllSlices2D); m_BtnApply3D = new QPushButton("Confirm", m_GroupBoxEnableExclusiveInterpolationMode); vboxLayout->addWidget(m_BtnApply3D); // T28261 // m_BtnSuggestPlane = new QPushButton("Suggest a plane", m_GroupBoxEnableExclusiveInterpolationMode); // vboxLayout->addWidget(m_BtnSuggestPlane); m_BtnReinit3DInterpolation = new QPushButton("Reinit Interpolation", m_GroupBoxEnableExclusiveInterpolationMode); vboxLayout->addWidget(m_BtnReinit3DInterpolation); m_ChkShowPositionNodes = new QCheckBox("Show Position Nodes", m_GroupBoxEnableExclusiveInterpolationMode); vboxLayout->addWidget(m_ChkShowPositionNodes); this->HideAllInterpolationControls(); connect(m_CmbInterpolation, SIGNAL(currentIndexChanged(int)), this, SLOT(OnInterpolationMethodChanged(int))); connect(m_BtnApply2D, SIGNAL(clicked()), this, SLOT(OnAcceptInterpolationClicked())); connect(m_BtnApplyForAllSlices2D, SIGNAL(clicked()), this, SLOT(OnAcceptAllInterpolationsClicked())); connect(m_BtnApply3D, SIGNAL(clicked()), this, SLOT(OnAccept3DInterpolationClicked())); connect(m_BtnReinit3DInterpolation, SIGNAL(clicked()), this, SLOT(OnReinit3DInterpolation())); connect(m_ChkShowPositionNodes, SIGNAL(toggled(bool)), this, SLOT(OnShowMarkers(bool))); connect(m_ChkShowPositionNodes, SIGNAL(toggled(bool)), this, SIGNAL(SignalShowMarkerNodes(bool))); QHBoxLayout *layout = new QHBoxLayout(this); layout->addWidget(m_GroupBoxEnableExclusiveInterpolationMode); this->setLayout(layout); itk::ReceptorMemberCommand::Pointer command = itk::ReceptorMemberCommand::New(); command->SetCallbackFunction(this, &QmitkSlicesInterpolator::OnInterpolationInfoChanged); InterpolationInfoChangedObserverTag = m_Interpolator->AddObserver(itk::ModifiedEvent(), command); itk::ReceptorMemberCommand::Pointer command2 = itk::ReceptorMemberCommand::New(); command2->SetCallbackFunction(this, &QmitkSlicesInterpolator::OnSurfaceInterpolationInfoChanged); SurfaceInterpolationInfoChangedObserverTag = m_SurfaceInterpolator->AddObserver(itk::ModifiedEvent(), command2); auto command3 = itk::ReceptorMemberCommand::New(); command3->SetCallbackFunction(this, &QmitkSlicesInterpolator::OnInterpolationAborted); InterpolationAbortedObserverTag = m_Interpolator->AddObserver(itk::AbortEvent(), command3); // feedback node and its visualization properties m_FeedbackNode = mitk::DataNode::New(); mitk::CoreObjectFactory::GetInstance()->SetDefaultProperties(m_FeedbackNode); m_FeedbackNode->SetProperty("binary", mitk::BoolProperty::New(true)); m_FeedbackNode->SetProperty("outline binary", mitk::BoolProperty::New(true)); m_FeedbackNode->SetProperty("color", mitk::ColorProperty::New(255.0, 255.0, 0.0)); m_FeedbackNode->SetProperty("texture interpolation", mitk::BoolProperty::New(false)); m_FeedbackNode->SetProperty("layer", mitk::IntProperty::New(20)); m_FeedbackNode->SetProperty("levelwindow", mitk::LevelWindowProperty::New(mitk::LevelWindow(0, 1))); m_FeedbackNode->SetProperty("name", mitk::StringProperty::New("Interpolation feedback")); m_FeedbackNode->SetProperty("opacity", mitk::FloatProperty::New(0.8)); m_FeedbackNode->SetProperty("helper object", mitk::BoolProperty::New(true)); m_InterpolatedSurfaceNode = mitk::DataNode::New(); m_InterpolatedSurfaceNode->SetProperty("color", mitk::ColorProperty::New(SURFACE_COLOR_RGB)); m_InterpolatedSurfaceNode->SetProperty("name", mitk::StringProperty::New("Surface Interpolation feedback")); m_InterpolatedSurfaceNode->SetProperty("opacity", mitk::FloatProperty::New(0.5)); m_InterpolatedSurfaceNode->SetProperty("line width", mitk::FloatProperty::New(4.0f)); m_InterpolatedSurfaceNode->SetProperty("includeInBoundingBox", mitk::BoolProperty::New(false)); m_InterpolatedSurfaceNode->SetProperty("helper object", mitk::BoolProperty::New(true)); m_InterpolatedSurfaceNode->SetVisibility(false); m_3DContourNode = mitk::DataNode::New(); m_3DContourNode->SetProperty("color", mitk::ColorProperty::New(0.0, 0.0, 0.0)); m_3DContourNode->SetProperty("hidden object", mitk::BoolProperty::New(true)); m_3DContourNode->SetProperty("name", mitk::StringProperty::New("Drawn Contours")); m_3DContourNode->SetProperty("material.representation", mitk::VtkRepresentationProperty::New(VTK_WIREFRAME)); m_3DContourNode->SetProperty("material.wireframeLineWidth", mitk::FloatProperty::New(2.0f)); m_3DContourNode->SetProperty("3DContourContainer", mitk::BoolProperty::New(true)); m_3DContourNode->SetProperty("includeInBoundingBox", mitk::BoolProperty::New(false)); m_3DContourNode->SetVisibility(false); QWidget::setContentsMargins(0, 0, 0, 0); if (QWidget::layout() != nullptr) { QWidget::layout()->setContentsMargins(0, 0, 0, 0); } // For running 3D Interpolation in background // create a QFuture and a QFutureWatcher connect(&m_Watcher, SIGNAL(started()), this, SLOT(StartUpdateInterpolationTimer())); connect(&m_Watcher, SIGNAL(finished()), this, SLOT(OnSurfaceInterpolationFinished())); connect(&m_Watcher, SIGNAL(finished()), this, SLOT(StopUpdateInterpolationTimer())); m_Timer = new QTimer(this); connect(m_Timer, SIGNAL(timeout()), this, SLOT(ChangeSurfaceColor())); } void QmitkSlicesInterpolator::SetDataStorage(mitk::DataStorage::Pointer storage) { if (m_DataStorage == storage) { return; } if (m_DataStorage.IsNotNull()) { m_DataStorage->RemoveNodeEvent.RemoveListener( mitk::MessageDelegate1(this, &QmitkSlicesInterpolator::NodeRemoved) ); } m_DataStorage = storage; m_SurfaceInterpolator->SetDataStorage(storage); if (m_DataStorage.IsNotNull()) { m_DataStorage->RemoveNodeEvent.AddListener( mitk::MessageDelegate1(this, &QmitkSlicesInterpolator::NodeRemoved) ); } } mitk::DataStorage *QmitkSlicesInterpolator::GetDataStorage() { if (m_DataStorage.IsNotNull()) { return m_DataStorage; } else { return nullptr; } } void QmitkSlicesInterpolator::InitializeWindow(QmitkRenderWindow* window) { auto slicer = window->GetSliceNavigationController(); if (slicer == nullptr) { MITK_WARN << "Tried setting up interpolation for a render window that does not have a slice navigation controller set"; return; } // Has to be initialized m_LastSNC = slicer; m_TimePoints.insert(slicer, slicer->GetSelectedTimePoint()); itk::MemberCommand::Pointer deleteCommand = itk::MemberCommand::New(); deleteCommand->SetCallbackFunction(this, &QmitkSlicesInterpolator::OnSliceNavigationControllerDeleted); m_ControllerToDeleteObserverTag[slicer] = slicer->AddObserver(itk::DeleteEvent(), deleteCommand); itk::MemberCommand::Pointer timeChangedCommand = itk::MemberCommand::New(); timeChangedCommand->SetCallbackFunction(this, &QmitkSlicesInterpolator::OnTimeChanged); m_ControllerToTimeObserverTag[slicer] = slicer->AddObserver(mitk::SliceNavigationController::TimeGeometryEvent(nullptr, 0), timeChangedCommand); itk::MemberCommand::Pointer sliceChangedCommand = itk::MemberCommand::New(); sliceChangedCommand->SetCallbackFunction(this, &QmitkSlicesInterpolator::OnSliceChanged); m_ControllerToSliceObserverTag[slicer] = slicer->AddObserver(mitk::SliceNavigationController::GeometrySliceEvent(nullptr, 0), sliceChangedCommand); } void QmitkSlicesInterpolator::Initialize(mitk::ToolManager *toolManager, const QList& windows) { Q_ASSERT(!windows.empty()); if (m_Initialized) { // remove old observers this->Uninitialize(); } m_ToolManager = toolManager; if (m_ToolManager) { // set enabled only if a segmentation is selected mitk::DataNode *node = m_ToolManager->GetWorkingData(0); QWidget::setEnabled(node != nullptr); // react whenever the set of selected segmentation changes m_ToolManager->WorkingDataChanged += mitk::MessageDelegate(this, &QmitkSlicesInterpolator::OnToolManagerWorkingDataModified); m_ToolManager->ReferenceDataChanged += mitk::MessageDelegate( this, &QmitkSlicesInterpolator::OnToolManagerReferenceDataModified); // connect to the slice navigation controller. after each change, call the interpolator for (auto* window : windows) { this->InitializeWindow(window); } m_ActionToSlicer = createActionToSlicer(windows); } m_Initialized = true; } void QmitkSlicesInterpolator::Uninitialize() { if (m_ToolManager.IsNotNull()) { m_ToolManager->WorkingDataChanged -= mitk::MessageDelegate(this, &QmitkSlicesInterpolator::OnToolManagerWorkingDataModified); m_ToolManager->ReferenceDataChanged -= mitk::MessageDelegate( this, &QmitkSlicesInterpolator::OnToolManagerReferenceDataModified); } for (auto* slicer : m_ControllerToTimeObserverTag.keys()) { slicer->RemoveObserver(m_ControllerToDeleteObserverTag.take(slicer)); slicer->RemoveObserver(m_ControllerToTimeObserverTag.take(slicer)); slicer->RemoveObserver(m_ControllerToSliceObserverTag.take(slicer)); } this->ClearSegmentationObservers(); m_ActionToSlicer.clear(); m_ToolManager = nullptr; m_Initialized = false; } QmitkSlicesInterpolator::~QmitkSlicesInterpolator() { if (m_Initialized) { // remove old observers this->Uninitialize(); } WaitForFutures(); if (m_DataStorage.IsNotNull()) { m_DataStorage->RemoveNodeEvent.RemoveListener( mitk::MessageDelegate1(this, &QmitkSlicesInterpolator::NodeRemoved) ); if (m_DataStorage->Exists(m_3DContourNode)) m_DataStorage->Remove(m_3DContourNode); if (m_DataStorage->Exists(m_InterpolatedSurfaceNode)) m_DataStorage->Remove(m_InterpolatedSurfaceNode); } // remove observer m_Interpolator->RemoveObserver(InterpolationAbortedObserverTag); m_Interpolator->RemoveObserver(InterpolationInfoChangedObserverTag); m_SurfaceInterpolator->RemoveObserver(SurfaceInterpolationInfoChangedObserverTag); m_SurfaceInterpolator->UnsetSelectedImage(); delete m_Timer; } /** External enableization... */ void QmitkSlicesInterpolator::setEnabled(bool enable) { QWidget::setEnabled(enable); // Set the gui elements of the different interpolation modi enabled if (enable) { if (m_2DInterpolationEnabled) { this->Show2DInterpolationControls(true); m_Interpolator->Activate2DInterpolation(true); } else if (m_3DInterpolationEnabled) { this->Show3DInterpolationControls(true); this->Show3DInterpolationResult(true); } } // Set all gui elements of the interpolation disabled else { this->HideAllInterpolationControls(); this->Show3DInterpolationResult(false); } } void QmitkSlicesInterpolator::On2DInterpolationEnabled(bool status) { OnInterpolationActivated(status); m_Interpolator->Activate2DInterpolation(status); } void QmitkSlicesInterpolator::On3DInterpolationEnabled(bool status) { On3DInterpolationActivated(status); } void QmitkSlicesInterpolator::OnInterpolationDisabled(bool status) { if (status) { OnInterpolationActivated(!status); On3DInterpolationActivated(!status); this->Show3DInterpolationResult(false); } } void QmitkSlicesInterpolator::HideAllInterpolationControls() { this->Show2DInterpolationControls(false); this->Show3DInterpolationControls(false); } void QmitkSlicesInterpolator::Show2DInterpolationControls(bool show) { m_BtnApply2D->setVisible(show); m_BtnApplyForAllSlices2D->setVisible(show); } void QmitkSlicesInterpolator::Show3DInterpolationControls(bool show) { m_BtnApply3D->setVisible(show); // T28261 // m_BtnSuggestPlane->setVisible(show); m_ChkShowPositionNodes->setVisible(show); m_BtnReinit3DInterpolation->setVisible(show); } void QmitkSlicesInterpolator::OnInterpolationMethodChanged(int index) { switch (index) { case 0: // Disabled m_GroupBoxEnableExclusiveInterpolationMode->setTitle("Interpolation"); this->HideAllInterpolationControls(); this->OnInterpolationActivated(false); this->On3DInterpolationActivated(false); this->Show3DInterpolationResult(false); m_Interpolator->Activate2DInterpolation(false); break; case 1: // 2D m_GroupBoxEnableExclusiveInterpolationMode->setTitle("Interpolation (Enabled)"); this->HideAllInterpolationControls(); this->Show2DInterpolationControls(true); this->OnInterpolationActivated(true); this->On3DInterpolationActivated(false); m_Interpolator->Activate2DInterpolation(true); break; case 2: // 3D m_GroupBoxEnableExclusiveInterpolationMode->setTitle("Interpolation (Enabled)"); this->HideAllInterpolationControls(); this->Show3DInterpolationControls(true); this->OnInterpolationActivated(false); this->On3DInterpolationActivated(true); m_Interpolator->Activate2DInterpolation(false); break; default: MITK_ERROR << "Unknown interpolation method!"; m_CmbInterpolation->setCurrentIndex(0); break; } } void QmitkSlicesInterpolator::OnShowMarkers(bool state) { mitk::DataStorage::SetOfObjects::ConstPointer allContourMarkers = m_DataStorage->GetSubset(mitk::NodePredicateProperty::New("isContourMarker", mitk::BoolProperty::New(true))); for (mitk::DataStorage::SetOfObjects::ConstIterator it = allContourMarkers->Begin(); it != allContourMarkers->End(); ++it) { it->Value()->SetProperty("helper object", mitk::BoolProperty::New(!state)); } } void QmitkSlicesInterpolator::OnToolManagerWorkingDataModified() { this->ClearSegmentationObservers(); if (m_ToolManager->GetWorkingData(0) != nullptr) { m_Segmentation = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); auto labelSetImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); m_BtnReinit3DInterpolation->setEnabled(true); try { if (m_SegmentationObserverTags.find(labelSetImage) == m_SegmentationObserverTags.end()) { auto command2 = itk::MemberCommand::New(); command2->SetCallbackFunction(this, &QmitkSlicesInterpolator::OnModifyLabelChanged); auto workingImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); m_SegmentationObserverTags[workingImage] = workingImage->AddObserver(itk::ModifiedEvent(), command2); } } catch (const std::exception& e) { MITK_ERROR << "Error casting node data to LabelSetImage\n"; } } else { // If no workingdata is set, remove the interpolation feedback this->GetDataStorage()->Remove(m_FeedbackNode); m_FeedbackNode->SetData(nullptr); this->GetDataStorage()->Remove(m_3DContourNode); m_3DContourNode->SetData(nullptr); this->GetDataStorage()->Remove(m_InterpolatedSurfaceNode); m_InterpolatedSurfaceNode->SetData(nullptr); m_BtnReinit3DInterpolation->setEnabled(false); m_CmbInterpolation->setCurrentIndex(0); return; } // Updating the current selected segmentation for the 3D interpolation this->SetCurrentContourListID(); if (m_2DInterpolationEnabled) { OnInterpolationActivated(true); // re-initialize if needed } } void QmitkSlicesInterpolator::OnToolManagerReferenceDataModified() { } void QmitkSlicesInterpolator::OnTimeChanged(itk::Object *sender, const itk::EventObject &e) { // Check if we really have a GeometryTimeEvent if (!dynamic_cast(&e)) return; mitk::SliceNavigationController *slicer = dynamic_cast(sender); Q_ASSERT(slicer); const auto timePoint = slicer->GetSelectedTimePoint(); m_TimePoints[slicer] = timePoint; if (m_Watcher.isRunning()) m_Watcher.waitForFinished(); if (timePoint != m_SurfaceInterpolator->GetCurrentTimePoint()) { m_SurfaceInterpolator->SetCurrentTimePoint(timePoint); if (m_3DInterpolationEnabled) { m_3DContourNode->SetData(nullptr); m_InterpolatedSurfaceNode->SetData(nullptr); } m_SurfaceInterpolator->Modified(); } if (m_LastSNC == slicer) { slicer->SendSlice(); // will trigger a new interpolation } } void QmitkSlicesInterpolator::OnSliceChanged(itk::Object *sender, const itk::EventObject &e) { // Check whether we really have a GeometrySliceEvent if (!dynamic_cast(&e)) return; mitk::SliceNavigationController *slicer = dynamic_cast(sender); if(m_2DInterpolationEnabled) { this->On2DInterpolationEnabled(m_2DInterpolationEnabled); } if (TranslateAndInterpolateChangedSlice(e, slicer)) { slicer->GetRenderer()->RequestUpdate(); } } bool QmitkSlicesInterpolator::TranslateAndInterpolateChangedSlice(const itk::EventObject &e, mitk::SliceNavigationController *slicer) { if (!m_2DInterpolationEnabled) return false; try { const mitk::SliceNavigationController::GeometrySliceEvent &event = dynamic_cast(e); mitk::TimeGeometry *tsg = event.GetTimeGeometry(); if (tsg && m_TimePoints.contains(slicer) && tsg->IsValidTimePoint(m_TimePoints[slicer])) { mitk::SlicedGeometry3D *slicedGeometry = dynamic_cast(tsg->GetGeometryForTimePoint(m_TimePoints[slicer]).GetPointer()); if (slicedGeometry) { m_LastSNC = slicer; mitk::PlaneGeometry *plane = dynamic_cast(slicedGeometry->GetPlaneGeometry(event.GetPos())); if (plane) { Interpolate(plane, m_TimePoints[slicer], slicer); } return true; } } } catch (const std::bad_cast &) { return false; // so what } return false; } void QmitkSlicesInterpolator::OnLayerChanged() { auto* workingNode = m_ToolManager->GetWorkingData(0); if (workingNode != nullptr) { m_3DContourNode->SetData(nullptr); this->Show3DInterpolationResult(false); } if (m_3DInterpolationEnabled) { m_SurfaceInterpolator->Modified(); } if (m_2DInterpolationEnabled) { m_FeedbackNode->SetData(nullptr); this->OnInterpolationActivated(true); m_LastSNC->SendSlice(); } mitk::RenderingManager::GetInstance()->RequestUpdateAll(); this->UpdateVisibleSuggestion(); } void QmitkSlicesInterpolator::Interpolate(mitk::PlaneGeometry *plane, mitk::TimePointType timePoint, mitk::SliceNavigationController *slicer) { if (m_ToolManager) { mitk::DataNode *node = m_ToolManager->GetWorkingData(0); if (node) { m_Segmentation = dynamic_cast(node->GetData()); if (m_Segmentation) { if (!m_Segmentation->GetTimeGeometry()->IsValidTimePoint(timePoint)) { MITK_WARN << "Cannot interpolate segmentation. Passed time point is not within the time bounds of WorkingImage. Time point: " << timePoint; return; } const auto timeStep = m_Segmentation->GetTimeGeometry()->TimePointToTimeStep(timePoint); int clickedSliceDimension = -1; int clickedSliceIndex = -1; // calculate real slice position, i.e. slice of the image mitk::SegTool2D::DetermineAffectedImageSlice(m_Segmentation, plane, clickedSliceDimension, clickedSliceIndex); mitk::Image::Pointer interpolation = m_Interpolator->Interpolate(clickedSliceDimension, clickedSliceIndex, plane, timeStep); m_FeedbackNode->SetData(interpolation); // maybe just have a variable that stores the active label color. if (m_ToolManager) { auto* workingNode = m_ToolManager->GetWorkingData(0); if (workingNode != nullptr) { auto* activeLabel = dynamic_cast(workingNode->GetData())->GetActiveLabelSet()->GetActiveLabel(); if (nullptr != activeLabel) { auto activeColor = activeLabel->GetColor(); m_FeedbackNode->SetProperty("color", mitk::ColorProperty::New(activeColor)); } } } m_LastSNC = slicer; m_LastSliceIndex = clickedSliceIndex; } } } } void QmitkSlicesInterpolator::OnSurfaceInterpolationFinished() { mitk::Surface::Pointer interpolatedSurface = m_SurfaceInterpolator->GetInterpolationResult(); mitk::DataNode *workingNode = m_ToolManager->GetWorkingData(0); mitk::PlaneGeometry::Pointer slicingPlane = mitk::PlaneGeometry::New(); mitk::Vector3D slicingPlaneNormalVector; FillVector3D(slicingPlaneNormalVector,0.0,1.0,0.0); mitk::Point3D origin; FillVector3D(origin, 0.0, 0.0, 0.0); slicingPlane->InitializePlane(origin, slicingPlaneNormalVector); if (interpolatedSurface.IsNotNull() && workingNode) { m_BtnApply3D->setEnabled(true); // T28261 // m_BtnSuggestPlane->setEnabled(true); m_InterpolatedSurfaceNode->SetData(interpolatedSurface); m_3DContourNode->SetData(m_SurfaceInterpolator->GetContoursAsSurface()); this->Show3DInterpolationResult(true); if (!m_DataStorage->Exists(m_InterpolatedSurfaceNode)) { m_DataStorage->Add(m_InterpolatedSurfaceNode); } } else if (interpolatedSurface.IsNull()) { m_BtnApply3D->setEnabled(false); // T28261 // m_BtnSuggestPlane->setEnabled(false); if (m_DataStorage->Exists(m_InterpolatedSurfaceNode)) { this->Show3DInterpolationResult(false); } } m_BtnReinit3DInterpolation->setEnabled(true); for (auto* slicer : m_ControllerToTimeObserverTag.keys()) { slicer->GetRenderer()->RequestUpdate(); } m_SurfaceInterpolator->ReinitializeInterpolation(); } void QmitkSlicesInterpolator::OnAcceptInterpolationClicked() { auto* workingNode = m_ToolManager->GetWorkingData(0); auto* planeGeometry = m_LastSNC->GetCurrentPlaneGeometry(); auto* interpolatedPreview = dynamic_cast(m_FeedbackNode->GetData()); if (nullptr == workingNode || nullptr == interpolatedPreview) return; auto* segmentationImage = dynamic_cast(workingNode->GetData()); if (nullptr == segmentationImage) return; const auto timePoint = m_LastSNC->GetSelectedTimePoint(); if (!segmentationImage->GetTimeGeometry()->IsValidTimePoint(timePoint)) { MITK_WARN << "Cannot accept interpolation. Time point selected by SliceNavigationController is not within the time bounds of segmentation. Time point: " << timePoint; return; } const auto timeStep = segmentationImage->GetTimeGeometry()->TimePointToTimeStep(timePoint); auto interpolatedSlice = mitk::SegTool2D::GetAffectedImageSliceAs2DImage(planeGeometry, segmentationImage, timeStep)->Clone(); auto labelSet = segmentationImage->GetActiveLabelSet(); auto activeValue = labelSet->GetActiveLabel()->GetValue(); mitk::TransferLabelContentAtTimeStep( interpolatedPreview, interpolatedSlice, labelSet, timeStep, 0, mitk::LabelSetImage::UnlabeledValue, false, { {0, mitk::LabelSetImage::UnlabeledValue}, {1, activeValue} } ); mitk::SegTool2D::WriteBackSegmentationResult(workingNode, planeGeometry, interpolatedSlice, timeStep); m_FeedbackNode->SetData(nullptr); } void QmitkSlicesInterpolator::AcceptAllInterpolations(mitk::SliceNavigationController *slicer) { /* * What exactly is done here: * 1. We create an empty diff image for the current segmentation * 2. All interpolated slices are written into the diff image * 3. Then the diffimage is applied to the original segmentation */ if (m_Segmentation) { mitk::Image::Pointer segmentation3D = m_Segmentation; unsigned int timeStep = 0; const auto timePoint = slicer->GetSelectedTimePoint(); if (4 == m_Segmentation->GetDimension()) { const auto* geometry = m_Segmentation->GetTimeGeometry(); if (!geometry->IsValidTimePoint(timePoint)) { MITK_WARN << "Cannot accept all interpolations. Time point selected by passed SliceNavigationController is not within the time bounds of segmentation. Time point: " << timePoint; return; } mitk::Image::Pointer activeLabelImage; try { auto labelSetImage = dynamic_cast(m_Segmentation); activeLabelImage = labelSetImage->CreateLabelMask(labelSetImage->GetActiveLabelSet()->GetActiveLabel()->GetValue(), true, 0); } catch (const std::exception& e) { MITK_ERROR << e.what() << " | NO LABELSETIMAGE IN WORKING NODE\n"; } m_Interpolator->SetSegmentationVolume(activeLabelImage); timeStep = geometry->TimePointToTimeStep(timePoint); auto timeSelector = mitk::ImageTimeSelector::New(); timeSelector->SetInput(m_Segmentation); timeSelector->SetTimeNr(timeStep); timeSelector->Update(); segmentation3D = timeSelector->GetOutput(); } // Create an empty diff image for the undo operation auto diffImage = mitk::Image::New(); diffImage->Initialize(segmentation3D); // Create scope for ImageWriteAccessor so that the accessor is destroyed right after use { mitk::ImageWriteAccessor accessor(diffImage); // Set all pixels to zero auto pixelType = mitk::MakeScalarPixelType(); // For legacy purpose support former pixel type of segmentations (before multilabel) if (itk::IOComponentEnum::UCHAR == m_Segmentation->GetImageDescriptor()->GetChannelDescriptor().GetPixelType().GetComponentType()) pixelType = mitk::MakeScalarPixelType(); memset(accessor.GetData(), 0, pixelType.GetSize() * diffImage->GetDimension(0) * diffImage->GetDimension(1) * diffImage->GetDimension(2)); } // Since we need to shift the plane it must be clone so that the original plane isn't altered auto slicedGeometry = m_Segmentation->GetSlicedGeometry(); auto planeGeometry = slicer->GetCurrentPlaneGeometry()->Clone(); int sliceDimension = -1; int sliceIndex = -1; mitk::SegTool2D::DetermineAffectedImageSlice(m_Segmentation, planeGeometry, sliceDimension, sliceIndex); const auto numSlices = m_Segmentation->GetDimension(sliceDimension); mitk::ProgressBar::GetInstance()->AddStepsToDo(numSlices); std::atomic_uint totalChangedSlices; // Reuse interpolation algorithm instance for each slice to cache boundary calculations auto algorithm = mitk::ShapeBasedInterpolationAlgorithm::New(); // Distribute slice interpolations to multiple threads const auto numThreads = std::min(std::thread::hardware_concurrency(), numSlices); // const auto numThreads = 1; std::vector> sliceIndices(numThreads); for (std::remove_const_t sliceIndex = 0; sliceIndex < numSlices; ++sliceIndex) sliceIndices[sliceIndex % numThreads].push_back(sliceIndex); std::vector threads; threads.reserve(numThreads); // This lambda will be executed by the threads auto interpolate = [=, &interpolator = m_Interpolator, &totalChangedSlices](unsigned int threadIndex) { auto clonedPlaneGeometry = planeGeometry->Clone(); auto origin = clonedPlaneGeometry->GetOrigin(); // Go through the sliced indices for (auto sliceIndex : sliceIndices[threadIndex]) { slicedGeometry->WorldToIndex(origin, origin); origin[sliceDimension] = sliceIndex; slicedGeometry->IndexToWorld(origin, origin); clonedPlaneGeometry->SetOrigin(origin); auto interpolation = interpolator->Interpolate(sliceDimension, sliceIndex, clonedPlaneGeometry, timeStep, algorithm); if (interpolation.IsNotNull()) { // Setting up the reslicing pipeline which allows us to write the interpolation results back into the image volume auto reslicer = vtkSmartPointer::New(); // Set overwrite mode to true to write back to the image volume reslicer->SetInputSlice(interpolation->GetSliceData()->GetVtkImageAccessor(interpolation)->GetVtkImageData()); reslicer->SetOverwriteMode(true); reslicer->Modified(); auto diffSliceWriter = mitk::ExtractSliceFilter::New(reslicer); diffSliceWriter->SetInput(diffImage); diffSliceWriter->SetTimeStep(0); diffSliceWriter->SetWorldGeometry(clonedPlaneGeometry); diffSliceWriter->SetVtkOutputRequest(true); diffSliceWriter->SetResliceTransformByGeometry(diffImage->GetTimeGeometry()->GetGeometryForTimeStep(0)); diffSliceWriter->Modified(); diffSliceWriter->Update(); ++totalChangedSlices; } mitk::ProgressBar::GetInstance()->Progress(); } }; m_Interpolator->EnableSliceImageCache(); // Do the interpolation here. for (size_t threadIndex = 0; threadIndex < numThreads; ++threadIndex) { interpolate(threadIndex); } m_Interpolator->DisableSliceImageCache(); const mitk::Label::PixelType newDestinationLabel = dynamic_cast(m_Segmentation)->GetActiveLabelSet()->GetActiveLabel()->GetValue(); // Do and Undo Operations if (totalChangedSlices > 0) { // Create do/undo operations auto* doOp = new mitk::ApplyDiffImageOperation(mitk::OpTEST, m_Segmentation, diffImage, timeStep); auto* undoOp = new mitk::ApplyDiffImageOperation(mitk::OpTEST, m_Segmentation, diffImage, timeStep); undoOp->SetFactor(-1.0); auto comment = "Confirm all interpolations (" + std::to_string(totalChangedSlices) + ")"; auto* undoStackItem = new mitk::OperationEvent(mitk::DiffImageApplier::GetInstanceForUndo(), doOp, undoOp, comment); mitk::OperationEvent::IncCurrGroupEventId(); mitk::OperationEvent::IncCurrObjectEventId(); mitk::UndoController::GetCurrentUndoModel()->SetOperationEvent(undoStackItem); mitk::DiffImageApplier::GetInstanceForUndo()->SetDestinationLabel(newDestinationLabel); // Apply the changes to the original image mitk::DiffImageApplier::GetInstanceForUndo()->ExecuteOperation(doOp); } m_FeedbackNode->SetData(nullptr); } mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkSlicesInterpolator::FinishInterpolation(mitk::SliceNavigationController *slicer) { // this redirect is for calling from outside if (slicer == nullptr) OnAcceptAllInterpolationsClicked(); else AcceptAllInterpolations(slicer); } void QmitkSlicesInterpolator::OnAcceptAllInterpolationsClicked() { QMenu orientationPopup(this); for (auto it = m_ActionToSlicer.begin(); it != m_ActionToSlicer.end(); ++it) orientationPopup.addAction(it->first); connect(&orientationPopup, SIGNAL(triggered(QAction *)), this, SLOT(OnAcceptAllPopupActivated(QAction *))); orientationPopup.exec(QCursor::pos()); } void QmitkSlicesInterpolator::OnAccept3DInterpolationClicked() { auto referenceImage = GetData(m_ToolManager->GetReferenceData(0)); auto* segmentationDataNode = m_ToolManager->GetWorkingData(0); auto labelSetImage = dynamic_cast(segmentationDataNode->GetData()); auto activeLabelColor = labelSetImage->GetActiveLabelSet()->GetActiveLabel()->GetColor(); std::string activeLabelName = labelSetImage->GetActiveLabelSet()->GetActiveLabel()->GetName(); auto segmentation = GetData(segmentationDataNode); if (referenceImage.IsNull() || segmentation.IsNull()) return; const auto* segmentationGeometry = segmentation->GetTimeGeometry(); const auto timePoint = m_LastSNC->GetSelectedTimePoint(); if (!referenceImage->GetTimeGeometry()->IsValidTimePoint(timePoint) || !segmentationGeometry->IsValidTimePoint(timePoint)) { MITK_WARN << "Cannot accept interpolation. Current time point is not within the time bounds of the patient image and segmentation."; return; } auto interpolatedSurface = GetData(m_InterpolatedSurfaceNode); if (interpolatedSurface.IsNull()) return; auto surfaceToImageFilter = mitk::SurfaceToImageFilter::New(); surfaceToImageFilter->SetImage(referenceImage); surfaceToImageFilter->SetMakeOutputBinary(true); surfaceToImageFilter->SetUShortBinaryPixelType(itk::IOComponentEnum::USHORT == segmentation->GetPixelType().GetComponentType()); surfaceToImageFilter->SetInput(interpolatedSurface); surfaceToImageFilter->Update(); mitk::Image::Pointer interpolatedSegmentation = surfaceToImageFilter->GetOutput(); auto timeStep = segmentationGeometry->TimePointToTimeStep(timePoint); const mitk::Label::PixelType newDestinationLabel = labelSetImage->GetActiveLabelSet()->GetActiveLabel()->GetValue(); TransferLabelContentAtTimeStep( interpolatedSegmentation, labelSetImage, labelSetImage->GetActiveLabelSet(), timeStep, 0, 0, false, {{1, newDestinationLabel}}, mitk::MultiLabelSegmentation::MergeStyle::Merge, mitk::MultiLabelSegmentation::OverwriteStyle::RegardLocks); // m_CmbInterpolation->setCurrentIndex(0); this->Show3DInterpolationResult(false); std::string name = segmentationDataNode->GetName() + " 3D-interpolation - " + activeLabelName; mitk::TimeBounds timeBounds; if (1 < interpolatedSurface->GetTimeSteps()) { name += "_t" + std::to_string(timeStep); auto* polyData = vtkPolyData::New(); polyData->DeepCopy(interpolatedSurface->GetVtkPolyData(timeStep)); auto surface = mitk::Surface::New(); surface->SetVtkPolyData(polyData); interpolatedSurface = surface; timeBounds = segmentationGeometry->GetTimeBounds(timeStep); } else { timeBounds = segmentationGeometry->GetTimeBounds(0); } auto* surfaceGeometry = static_cast(interpolatedSurface->GetTimeGeometry()); surfaceGeometry->SetFirstTimePoint(timeBounds[0]); surfaceGeometry->SetStepDuration(timeBounds[1] - timeBounds[0]); // Typical file formats for surfaces do not save any time-related information. As a workaround at least for MITK scene files, we have the // possibility to seralize this information as properties. interpolatedSurface->SetProperty("ProportionalTimeGeometry.FirstTimePoint", mitk::FloatProperty::New(surfaceGeometry->GetFirstTimePoint())); interpolatedSurface->SetProperty("ProportionalTimeGeometry.StepDuration", mitk::FloatProperty::New(surfaceGeometry->GetStepDuration())); auto interpolatedSurfaceDataNode = mitk::DataNode::New(); interpolatedSurfaceDataNode->SetData(interpolatedSurface); interpolatedSurfaceDataNode->SetName(name); interpolatedSurfaceDataNode->SetOpacity(0.7f); interpolatedSurfaceDataNode->SetColor(activeLabelColor); m_DataStorage->Add(interpolatedSurfaceDataNode, segmentationDataNode); } void QmitkSlicesInterpolator::OnReinit3DInterpolation() { // Step 1. Load from the isContourPlaneGeometry nodes the contourNodes. mitk::NodePredicateProperty::Pointer pred = mitk::NodePredicateProperty::New("isContourPlaneGeometry", mitk::BoolProperty::New(true)); mitk::DataStorage::SetOfObjects::ConstPointer contourNodes = m_DataStorage->GetDerivations(m_ToolManager->GetWorkingData(0), pred); if (contourNodes->Size() != 0) { std::vector contourPlanes; std::vector contourList; if (m_ToolManager->GetWorkingData(0) != nullptr) { try { auto labelSetImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); auto activeLayerID = labelSetImage->GetActiveLayer(); const auto timePoint = m_LastSNC->GetSelectedTimePoint(); if (!labelSetImage->GetTimeGeometry()->IsValidTimePoint(timePoint)) { MITK_ERROR << "Invalid time point requested for interpolation pipeline."; return; } // Adding layer, label and timeStep information for the contourNodes. for (auto it = contourNodes->Begin(); it != contourNodes->End(); ++it) { auto contourNode = it->Value(); auto layerID = dynamic_cast(contourNode->GetProperty("layerID"))->GetValue(); auto labelID = dynamic_cast(contourNode->GetProperty("labelID"))->GetValue(); auto timeStep = dynamic_cast(contourNode->GetProperty("timeStep"))->GetValue(); auto px = dynamic_cast(contourNode->GetProperty("px"))->GetValue(); auto py = dynamic_cast(contourNode->GetProperty("py"))->GetValue(); auto pz = dynamic_cast(contourNode->GetProperty("pz"))->GetValue(); // auto layerImage = labelSetImage->GetLayerImage(layerID); auto planeGeometry = dynamic_cast(contourNode->GetData())->GetPlaneGeometry(); labelSetImage->SetActiveLayer(layerID); auto sliceImage = ExtractSliceFromImage(labelSetImage, planeGeometry, timeStep); labelSetImage->SetActiveLayer(activeLayerID); mitk::ImageToContourFilter::Pointer contourExtractor = mitk::ImageToContourFilter::New(); contourExtractor->SetInput(sliceImage); contourExtractor->SetContourValue(labelID); contourExtractor->Update(); mitk::Surface::Pointer contour = contourExtractor->GetOutput(); if (contour->GetVtkPolyData()->GetNumberOfPoints() == 0) continue; vtkSmartPointer intArray = vtkSmartPointer::New(); intArray->InsertNextValue(labelID); intArray->InsertNextValue(layerID); intArray->InsertNextValue(timeStep); contour->GetVtkPolyData()->GetFieldData()->AddArray(intArray); vtkSmartPointer doubleArray = vtkSmartPointer::New(); doubleArray->InsertNextValue(px); doubleArray->InsertNextValue(py); doubleArray->InsertNextValue(pz); contour->GetVtkPolyData()->GetFieldData()->AddArray(doubleArray); contour->DisconnectPipeline(); contourList.push_back(contour); contourPlanes.push_back(planeGeometry); } labelSetImage->SetActiveLayer(activeLayerID); // size_t activeLayer = labelSetImage->GetActiveLayer(); for (size_t l = 0; l < labelSetImage->GetNumberOfLayers(); ++l) { this->OnAddLabelSetConnection(l); } // labelSetImage->SetActiveLayer(activeLayer); m_SurfaceInterpolator->CompleteReinitialization(contourList, contourPlanes); } catch(const std::exception& e) { MITK_ERROR << "Exception thrown casting toolmanager working data to labelsetImage"; } } } else { m_BtnApply3D->setEnabled(false); QMessageBox errorInfo; errorInfo.setWindowTitle("Reinitialize surface interpolation"); errorInfo.setIcon(QMessageBox::Information); errorInfo.setText("No contours available for the selected segmentation!"); errorInfo.exec(); } } void QmitkSlicesInterpolator::OnAcceptAllPopupActivated(QAction *action) { try { auto iter = m_ActionToSlicer.find(action); if (iter != m_ActionToSlicer.end()) { mitk::SliceNavigationController *slicer = iter->second; AcceptAllInterpolations(slicer); } } catch (...) { /* Showing message box with possible memory error */ QMessageBox errorInfo; errorInfo.setWindowTitle("Interpolation Process"); errorInfo.setIcon(QMessageBox::Critical); errorInfo.setText("An error occurred during interpolation. Possible cause: Not enough memory!"); errorInfo.exec(); // additional error message on std::cerr std::cerr << "Ill construction in " __FILE__ " l. " << __LINE__ << std::endl; } } void QmitkSlicesInterpolator::OnInterpolationActivated(bool on) { m_2DInterpolationEnabled = on; try { if (m_DataStorage.IsNotNull()) { if (on && !m_DataStorage->Exists(m_FeedbackNode)) { m_DataStorage->Add(m_FeedbackNode); } } } catch (...) { // don't care (double add/remove) } if (m_ToolManager) { mitk::DataNode *workingNode = m_ToolManager->GetWorkingData(0); mitk::DataNode *referenceNode = m_ToolManager->GetReferenceData(0); QWidget::setEnabled(workingNode != nullptr); m_BtnApply2D->setEnabled(on); m_FeedbackNode->SetVisibility(on); if (!on) { mitk::RenderingManager::GetInstance()->RequestUpdateAll(); return; } if (workingNode) { auto labelSetImage = dynamic_cast(workingNode->GetData()); if (nullptr == labelSetImage) { MITK_ERROR << "NO LABELSETIMAGE IN WORKING NODE\n"; mitk::RenderingManager::GetInstance()->RequestUpdateAll(); return; } auto* activeLabel = labelSetImage->GetActiveLabelSet()->GetActiveLabel(); auto* segmentation = dynamic_cast(workingNode->GetData()); if (nullptr != activeLabel && nullptr != segmentation) { auto activeLabelImage = labelSetImage->CreateLabelMask(activeLabel->GetValue(), true, 0); m_Interpolator->SetSegmentationVolume(activeLabelImage); if (referenceNode) { mitk::Image *referenceImage = dynamic_cast(referenceNode->GetData()); m_Interpolator->SetReferenceVolume(referenceImage); // may be nullptr } } } } this->UpdateVisibleSuggestion(); } void QmitkSlicesInterpolator::Run3DInterpolation() { m_SurfaceInterpolator->Interpolate(); } void QmitkSlicesInterpolator::StartUpdateInterpolationTimer() { m_Timer->start(500); } void QmitkSlicesInterpolator::StopUpdateInterpolationTimer() { if(m_ToolManager) { auto* workingNode = m_ToolManager->GetWorkingData(0); auto activeColor = dynamic_cast(workingNode->GetData())->GetActiveLabelSet()->GetActiveLabel()->GetColor(); m_InterpolatedSurfaceNode->SetProperty("color", mitk::ColorProperty::New(activeColor)); m_3DContourNode->SetProperty("color", mitk::ColorProperty::New(activeColor)); } m_Timer->stop(); } void QmitkSlicesInterpolator::ChangeSurfaceColor() { float currentColor[3]; m_InterpolatedSurfaceNode->GetColor(currentColor); m_InterpolatedSurfaceNode->SetProperty("color", mitk::ColorProperty::New(SURFACE_COLOR_RGB)); m_InterpolatedSurfaceNode->Update(); mitk::RenderingManager::GetInstance()->RequestUpdateAll(mitk::RenderingManager::REQUEST_UPDATE_3DWINDOWS); } void QmitkSlicesInterpolator::PrepareInputsFor3DInterpolation() { if (m_DataStorage.IsNotNull() && m_ToolManager && m_3DInterpolationEnabled) { auto *workingNode = m_ToolManager->GetWorkingData(0); if (workingNode != nullptr) { int ret = QMessageBox::Yes; if (m_SurfaceInterpolator->EstimatePortionOfNeededMemory() > 0.5) { QMessageBox msgBox; msgBox.setText("Due to short handed system memory the 3D interpolation may be very slow!"); msgBox.setInformativeText("Are you sure you want to activate the 3D interpolation?"); msgBox.setStandardButtons(QMessageBox::No | QMessageBox::Yes); ret = msgBox.exec(); } auto labelSetImage = dynamic_cast(workingNode->GetData()); auto activeLabel = labelSetImage->GetActiveLabelSet()->GetActiveLabel()->GetValue(); m_SurfaceInterpolator->AddActiveLabelContoursForInterpolation(activeLabel); if (m_Watcher.isRunning()) m_Watcher.waitForFinished(); if (ret == QMessageBox::Yes) { // Maybe set the segmentation node here m_Future = QtConcurrent::run(this, &QmitkSlicesInterpolator::Run3DInterpolation); m_Watcher.setFuture(m_Future); } else { m_CmbInterpolation->setCurrentIndex(0); } } else { QWidget::setEnabled(false); m_ChkShowPositionNodes->setEnabled(m_3DInterpolationEnabled); } } if (!m_3DInterpolationEnabled) { this->Show3DInterpolationResult(false); m_BtnApply3D->setEnabled(m_3DInterpolationEnabled); // T28261 // m_BtnSuggestPlane->setEnabled(m_3DInterpolationEnabled); } mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkSlicesInterpolator::On3DInterpolationActivated(bool on) { m_3DInterpolationEnabled = on; try { // this->PrepareInputsFor3DInterpolation(); m_SurfaceInterpolator->Modified(); } catch (...) { MITK_ERROR << "Error with 3D surface interpolation!"; } mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkSlicesInterpolator::EnableInterpolation(bool on) { // only to be called from the outside world // just a redirection to OnInterpolationActivated OnInterpolationActivated(on); } void QmitkSlicesInterpolator::Enable3DInterpolation(bool on) { // only to be called from the outside world // just a redirection to OnInterpolationActivated this->On3DInterpolationActivated(on); } void QmitkSlicesInterpolator::UpdateVisibleSuggestion() { mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkSlicesInterpolator::OnInterpolationInfoChanged(const itk::EventObject & /*e*/) { // something (e.g. undo) changed the interpolation info, we should refresh our display this->UpdateVisibleSuggestion(); } void QmitkSlicesInterpolator::OnInterpolationAborted(const itk::EventObject& /*e*/) { m_CmbInterpolation->setCurrentIndex(0); m_FeedbackNode->SetData(nullptr); } void QmitkSlicesInterpolator::OnSurfaceInterpolationInfoChanged(const itk::EventObject & /*e*/) { if (m_Watcher.isRunning()) m_Watcher.waitForFinished(); if (m_3DInterpolationEnabled) { m_3DContourNode->SetData(nullptr); m_InterpolatedSurfaceNode->SetData(nullptr); auto* workingNode = m_ToolManager->GetWorkingData(0); if (workingNode == nullptr) return; auto* labelSetImage = dynamic_cast(workingNode->GetData()); auto* label = labelSetImage->GetActiveLabelSet()->GetActiveLabel(); if (label == nullptr) return; m_SurfaceInterpolator->AddActiveLabelContoursForInterpolation(label->GetValue()); m_Future = QtConcurrent::run(this, &QmitkSlicesInterpolator::Run3DInterpolation); m_Watcher.setFuture(m_Future); } } void QmitkSlicesInterpolator::SetCurrentContourListID() { // New ContourList = hide current interpolation Show3DInterpolationResult(false); if (m_DataStorage.IsNotNull() && m_ToolManager && m_LastSNC) { mitk::DataNode *workingNode = m_ToolManager->GetWorkingData(0); try{ auto labelSetImage = dynamic_cast(workingNode->GetData()); for (size_t layerID = 0; layerID < labelSetImage->GetNumberOfLayers(); ++layerID) { this->OnAddLabelSetConnection(layerID); } } catch (std::exception &e) { MITK_ERROR << e.what() << "\n"; } if (workingNode) { QWidget::setEnabled(true); const auto timePoint = m_LastSNC->GetSelectedTimePoint(); // In case the time is not valid use 0 to access the time geometry of the working node unsigned int time_position = 0; if (!workingNode->GetData()->GetTimeGeometry()->IsValidTimePoint(timePoint)) { MITK_WARN << "Cannot accept interpolation. Time point selected by SliceNavigationController is not within the time bounds of WorkingImage. Time point: " << timePoint; return; } // Sets up the surface interpolator to accept time_position = workingNode->GetData()->GetTimeGeometry()->TimePointToTimeStep(timePoint); mitk::Vector3D spacing = workingNode->GetData()->GetGeometry(time_position)->GetSpacing(); double minSpacing = 100; double maxSpacing = 0; for (int i = 0; i < 3; i++) { if (spacing[i] < minSpacing) { minSpacing = spacing[i]; } if (spacing[i] > maxSpacing) { maxSpacing = spacing[i]; } } m_SurfaceInterpolator->SetMaxSpacing(maxSpacing); m_SurfaceInterpolator->SetMinSpacing(minSpacing); m_SurfaceInterpolator->SetDistanceImageVolume(50000); mitk::Image::Pointer segmentationImage; segmentationImage = dynamic_cast(workingNode->GetData()); m_SurfaceInterpolator->SetCurrentInterpolationSession(segmentationImage); m_SurfaceInterpolator->SetCurrentTimePoint(timePoint); } else { QWidget::setEnabled(false); } } } void QmitkSlicesInterpolator::Show3DInterpolationResult(bool status) { if (m_InterpolatedSurfaceNode.IsNotNull()) m_InterpolatedSurfaceNode->SetVisibility(status); if (m_3DContourNode.IsNotNull()) { auto allRenderWindows = mitk::BaseRenderer::GetAll3DRenderWindows(); for (auto mapit = allRenderWindows.begin(); mapit != allRenderWindows.end(); ++mapit) { m_3DContourNode->SetVisibility(status, mapit->second); } } mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkSlicesInterpolator::OnActiveLabelChanged(mitk::Label::PixelType) { m_3DContourNode->SetData(nullptr); m_FeedbackNode->SetData(nullptr); m_InterpolatedSurfaceNode->SetData(nullptr); if (m_Watcher.isRunning()) m_Watcher.waitForFinished(); if (m_3DInterpolationEnabled) { m_SurfaceInterpolator->Modified(); } if (m_2DInterpolationEnabled) { m_FeedbackNode->SetData(nullptr); this->OnInterpolationActivated(true); m_LastSNC->SendSlice(); } mitk::RenderingManager::GetInstance()->RequestUpdateAll(); this->UpdateVisibleSuggestion(); } void QmitkSlicesInterpolator::CheckSupportedImageDimension() { if (m_ToolManager->GetWorkingData(0)) { m_Segmentation = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); if (m_3DInterpolationEnabled && m_Segmentation && ((m_Segmentation->GetDimension() != 3) || (m_Segmentation->GetDimension() != 4)) ) { QMessageBox info; info.setWindowTitle("3D Interpolation Process"); info.setIcon(QMessageBox::Information); info.setText("3D Interpolation is only supported for 3D/4D images at the moment!"); info.exec(); m_CmbInterpolation->setCurrentIndex(0); } } } void QmitkSlicesInterpolator::OnSliceNavigationControllerDeleted(const itk::Object *sender, const itk::EventObject & /*e*/) { // Don't know how to avoid const_cast here?! mitk::SliceNavigationController *slicer = dynamic_cast(const_cast(sender)); if (slicer) { m_ControllerToTimeObserverTag.remove(slicer); m_ControllerToSliceObserverTag.remove(slicer); m_ControllerToDeleteObserverTag.remove(slicer); } } void QmitkSlicesInterpolator::WaitForFutures() { if (m_Watcher.isRunning()) { m_Watcher.waitForFinished(); } if (m_PlaneWatcher.isRunning()) { m_PlaneWatcher.waitForFinished(); } } void QmitkSlicesInterpolator::NodeRemoved(const mitk::DataNode* node) { if ((m_ToolManager && m_ToolManager->GetWorkingData(0) == node) || node == m_3DContourNode || node == m_FeedbackNode || node == m_InterpolatedSurfaceNode) { WaitForFutures(); } } void QmitkSlicesInterpolator::OnAddLabelSetConnection(unsigned int layerID) { if (m_ToolManager->GetWorkingData(0) != nullptr) { try { auto workingImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); auto labelSet = workingImage->GetLabelSet(layerID); labelSet->RemoveLabelEvent += mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnRemoveLabel); labelSet->ActiveLabelEvent += mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnActiveLabelChanged); workingImage->AfterChangeLayerEvent += mitk::MessageDelegate( this, &QmitkSlicesInterpolator::OnLayerChanged); m_SurfaceInterpolator->AddLabelSetConnection(layerID); } catch(const std::exception& e) { MITK_ERROR << e.what() << '\n'; } } } void QmitkSlicesInterpolator::OnAddLabelSetConnection() { if (m_ToolManager->GetWorkingData(0) != nullptr) { try { auto workingImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); workingImage->GetActiveLabelSet()->RemoveLabelEvent += mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnRemoveLabel); workingImage->GetActiveLabelSet()->ActiveLabelEvent += mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnActiveLabelChanged); workingImage->AfterChangeLayerEvent += mitk::MessageDelegate( this, &QmitkSlicesInterpolator::OnLayerChanged); m_SurfaceInterpolator->AddLabelSetConnection(); } catch(const std::exception& e) { MITK_ERROR << e.what() << '\n'; } } } void QmitkSlicesInterpolator::OnRemoveLabelSetConnection(mitk::LabelSetImage* labelSetImage, unsigned int layerID) { size_t previousLayerID = labelSetImage->GetActiveLayer(); labelSetImage->SetActiveLayer(layerID); labelSetImage->GetActiveLabelSet()->RemoveLabelEvent -= mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnRemoveLabel); labelSetImage->GetActiveLabelSet()->ActiveLabelEvent -= mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnActiveLabelChanged); labelSetImage->AfterChangeLayerEvent -= mitk::MessageDelegate( this, &QmitkSlicesInterpolator::OnLayerChanged); m_SurfaceInterpolator->RemoveLabelSetConnection(labelSetImage, layerID); labelSetImage->SetActiveLayer(previousLayerID); } void QmitkSlicesInterpolator::OnRemoveLabelSetConnection() { if (m_ToolManager->GetWorkingData(0) != nullptr) { try { auto workingImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); workingImage->GetActiveLabelSet()->RemoveLabelEvent -= mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnRemoveLabel); workingImage->GetActiveLabelSet()->ActiveLabelEvent -= mitk::MessageDelegate1( this, &QmitkSlicesInterpolator::OnActiveLabelChanged); workingImage->AfterChangeLayerEvent -= mitk::MessageDelegate( this, &QmitkSlicesInterpolator::OnLayerChanged); } catch(const std::exception& e) { MITK_ERROR << e.what() << '\n'; } } } void QmitkSlicesInterpolator::OnRemoveLabel(mitk::Label::PixelType /*removedLabelValue*/) { if (m_ToolManager->GetWorkingData(0) != nullptr) { try { auto labelSetImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); auto currentLayerID = labelSetImage->GetActiveLayer(); auto numTimeSteps = labelSetImage->GetTimeGeometry()->CountTimeSteps(); for (size_t t = 0; t < numTimeSteps; ++t) { m_SurfaceInterpolator->RemoveContours(m_PreviousActiveLabelValue,t,currentLayerID); } } catch(const std::exception& e) { MITK_ERROR << "Bad cast error for labelSetImage"; } } } void QmitkSlicesInterpolator::OnModifyLabelChanged(const itk::Object *caller, const itk::EventObject & /*event*/) { auto *tempImage = dynamic_cast(const_cast(caller) ) ; if( tempImage == nullptr) { MITK_ERROR << "Unable to cast caller to LabelSetImage."; return; } ModifyLabelActionTrigerred actionTriggered = ModifyLabelActionTrigerred::Null; if(m_ToolManager->GetWorkingData(0) != nullptr) { auto labelSetImage = dynamic_cast(m_ToolManager->GetWorkingData(0)->GetData()); if (labelSetImage == tempImage) { const auto timePoint = m_LastSNC->GetSelectedTimePoint(); if (!labelSetImage->GetTimeGeometry()->IsValidTimePoint(timePoint)) { MITK_ERROR << "Invalid time point requested for interpolation pipeline."; return; } auto timeStep = labelSetImage->GetTimeGeometry()->TimePointToTimeStep(timePoint); auto numLayersInCurrentSegmentation = m_SurfaceInterpolator->GetNumberOfLayersInCurrentSegmentation(); // This handles the add layer or remove layer operation. if (labelSetImage->GetNumberOfLayers() != numLayersInCurrentSegmentation) { bool addLayer = (labelSetImage->GetNumberOfLayers() == (numLayersInCurrentSegmentation +1) ); bool removeLayer = (labelSetImage->GetNumberOfLayers() == (numLayersInCurrentSegmentation - 1) ); m_SurfaceInterpolator->SetNumberOfLayersInCurrentSegmentation(labelSetImage->GetNumberOfLayers()); if (addLayer) { m_SurfaceInterpolator->OnAddLayer(); this->OnAddLabelSetConnection(); } if (removeLayer) { m_SurfaceInterpolator->OnRemoveLayer(); } return; } // Get the pixels present in the image. // This portion of the code deals with the merge and erase labels operations. auto imageDimension = labelSetImage->GetDimension(); if (imageDimension == 4) { actionTriggered = ModifyLabelProcessing<4>(labelSetImage, m_SurfaceInterpolator, timeStep); } else { actionTriggered = ModifyLabelProcessing<3>(labelSetImage, m_SurfaceInterpolator, timeStep); } if (actionTriggered == ModifyLabelActionTrigerred::Erase) { m_InterpolatedSurfaceNode->SetData(nullptr); } auto currentLayerID = labelSetImage->GetActiveLayer(); if (actionTriggered == ModifyLabelActionTrigerred::Merge) { this->MergeContours(timeStep, currentLayerID); m_SurfaceInterpolator->Modified(); } } } } void QmitkSlicesInterpolator::MergeContours(unsigned int timeStep, unsigned int layerID) { - std::vector& contours = - m_SurfaceInterpolator->GetContours(timeStep,layerID); + auto* contours = m_SurfaceInterpolator->GetContours(timeStep, layerID); + + if (nullptr == contours) + return; + std::this_thread::sleep_for(std::chrono::milliseconds(1000)); - for (size_t i = 0; i < contours.size(); ++i) + for (size_t i = 0; i < contours->size(); ++i) { - for (size_t j = i+1; j < contours.size(); ++j) + for (size_t j = i+1; j < contours->size(); ++j) { // And Labels are the same and Layers are the same. - bool areContoursCoplanar = AreContoursCoplanar(contours[i],contours[j]); + bool areContoursCoplanar = AreContoursCoplanar((*contours)[i], (*contours)[j]); - if ( areContoursCoplanar && (contours[i].LabelValue == contours[j].LabelValue) ) + if ( areContoursCoplanar && ((*contours)[i].LabelValue == (*contours)[j].LabelValue) ) { // Update the contour by re-extracting the slice from the corresponding plane. - mitk::Image::Pointer slice = ExtractSliceFromImage(m_Segmentation, contours[i].Plane, timeStep); + mitk::Image::Pointer slice = ExtractSliceFromImage(m_Segmentation, (*contours)[i].Plane, timeStep); mitk::ImageToContourFilter::Pointer contourExtractor = mitk::ImageToContourFilter::New(); contourExtractor->SetInput(slice); - contourExtractor->SetContourValue(contours[i].LabelValue); + contourExtractor->SetContourValue((*contours)[i].LabelValue); contourExtractor->Update(); mitk::Surface::Pointer contour = contourExtractor->GetOutput(); - contours[i].Contour = contour; + (*contours)[i].Contour = contour; // Update the interior point of the contour - contours[i].ContourPoint = m_SurfaceInterpolator->ComputeInteriorPointOfContour(contours[i],dynamic_cast(m_Segmentation)); + (*contours)[i].ContourPoint = m_SurfaceInterpolator->ComputeInteriorPointOfContour((*contours)[i],dynamic_cast(m_Segmentation)); // Setting the contour polygon data to an empty vtkPolyData, // as source label is empty after merge operation. - contours[j].Contour->SetVtkPolyData(vtkSmartPointer::New()); + (*contours)[j].Contour->SetVtkPolyData(vtkSmartPointer::New()); } } } auto segmentationNode = m_SurfaceInterpolator->GetSegmentationImageNode(); if (segmentationNode == nullptr) { MITK_ERROR << "segmentation Image Node not found\n"; } auto isContourPlaneGeometry = mitk::NodePredicateProperty::New("isContourPlaneGeometry", mitk::BoolProperty::New(true)); mitk::DataStorage::SetOfObjects::ConstPointer contourNodes = m_DataStorage->GetDerivations(segmentationNode, isContourPlaneGeometry); // Remove empty contour nodes. auto isContourEmpty = [] (const mitk::SurfaceInterpolationController::ContourPositionInformation& contour) { return (contour.Contour->GetVtkPolyData()->GetNumberOfPoints() == 0); }; - auto it = std::remove_if(contours.begin(), contours.end(), isContourEmpty); - contours.erase(it, contours.end()); + auto it = std::remove_if((*contours).begin(), (*contours).end(), isContourEmpty); + (*contours).erase(it, (*contours).end()); } void QmitkSlicesInterpolator::ClearSegmentationObservers() { auto dataIter = m_SegmentationObserverTags.begin(); while (dataIter != m_SegmentationObserverTags.end()) { auto labelSetImage = (*dataIter).first; labelSetImage->RemoveObserver((*dataIter).second); for (size_t layerID = 0; layerID < labelSetImage->GetNumberOfLayers(); ++layerID) { this->OnRemoveLabelSetConnection(labelSetImage, layerID); } ++dataIter; } m_SegmentationObserverTags.clear(); } diff --git a/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.cpp b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.cpp index 015595cda7..374f0ac5bd 100644 --- a/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.cpp +++ b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.cpp @@ -1,1402 +1,1399 @@ /*============================================================================ The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center (DKFZ) All rights reserved. Use of this source code is governed by a 3-clause BSD license that can be found in the LICENSE file. ============================================================================*/ #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include // Check whether the given contours are coplanar bool ContoursCoplanar(mitk::SurfaceInterpolationController::ContourPositionInformation leftHandSide, mitk::SurfaceInterpolationController::ContourPositionInformation rightHandSide) { // Here we check two things: // 1. Whether the normals of both contours are at least parallel // 2. Whether both contours lie in the same plane // Check for coplanarity: // a. Span a vector between two points one from each contour // b. Calculate dot product for the vector and one of the normals // c. If the dot is zero the two vectors are orthogonal and the contours are coplanar double vec[3]; vec[0] = leftHandSide.ContourPoint[0] - rightHandSide.ContourPoint[0]; vec[1] = leftHandSide.ContourPoint[1] - rightHandSide.ContourPoint[1]; vec[2] = leftHandSide.ContourPoint[2] - rightHandSide.ContourPoint[2]; double n[3]; n[0] = rightHandSide.ContourNormal[0]; n[1] = rightHandSide.ContourNormal[1]; n[2] = rightHandSide.ContourNormal[2]; double dot = vtkMath::Dot(n, vec); double n2[3]; n2[0] = leftHandSide.ContourNormal[0]; n2[1] = leftHandSide.ContourNormal[1]; n2[2] = leftHandSide.ContourNormal[2]; // The normals of both contours have to be parallel but not of the same orientation double lengthLHS = leftHandSide.ContourNormal.GetNorm(); double lengthRHS = rightHandSide.ContourNormal.GetNorm(); double dot2 = vtkMath::Dot(n, n2); bool contoursParallel = mitk::Equal(fabs(lengthLHS * lengthRHS), fabs(dot2), 0.001); if (mitk::Equal(dot, 0.0, 0.001) && contoursParallel) return true; else return false; } mitk::SurfaceInterpolationController::ContourPositionInformation CreateContourPositionInformation( mitk::Surface::Pointer contour, const mitk::PlaneGeometry* planeGeometry) { mitk::SurfaceInterpolationController::ContourPositionInformation contourInfo; contourInfo.Contour = contour; mitk::ScalarType n[3]; vtkPolygon::ComputeNormal(contour->GetVtkPolyData()->GetPoints(), n); contourInfo.ContourNormal = n; contourInfo.Pos = -1; contourInfo.TimeStep = std::numeric_limits::max(); contourInfo.Plane = const_cast(planeGeometry); auto contourIntArray = vtkIntArray::SafeDownCast( contour->GetVtkPolyData()->GetFieldData()->GetAbstractArray(0) ); if (contourIntArray->GetSize() < 2) { MITK_ERROR << "In CreateContourPositionInformation. The contourIntArray is empty."; } contourInfo.LabelValue = contourIntArray->GetValue(0); contourInfo.LayerValue = contourIntArray->GetValue(1); if (contourIntArray->GetSize() >= 3) { contourInfo.TimeStep = contourIntArray->GetValue(2); } contourInfo.SliceIndex = 0; return contourInfo; }; mitk::SurfaceInterpolationController::SurfaceInterpolationController() : m_SelectedSegmentation(nullptr), m_CurrentTimePoint(0.), m_ContourIndex(0), m_ContourPosIndex(0), m_NumberOfLayersInCurrentSegmentation(0), m_PreviousActiveLabelValue(0), m_CurrentActiveLabelValue(0), m_PreviousLayerIndex(0), m_CurrentLayerIndex(0) { m_DistanceImageSpacing = 0.0; m_ReduceFilter = ReduceContourSetFilter::New(); m_NormalsFilter = ComputeContourSetNormalsFilter::New(); m_InterpolateSurfaceFilter = CreateDistanceImageFromSurfaceFilter::New(); // m_TimeSelector = ImageTimeSelector::New(); m_ReduceFilter->SetUseProgressBar(false); // m_ReduceFilter->SetProgressStepSize(1); m_NormalsFilter->SetUseProgressBar(true); m_NormalsFilter->SetProgressStepSize(1); m_InterpolateSurfaceFilter->SetUseProgressBar(true); m_InterpolateSurfaceFilter->SetProgressStepSize(7); m_Contours = Surface::New(); m_PolyData = vtkSmartPointer::New(); vtkSmartPointer points = vtkSmartPointer::New(); m_PolyData->SetPoints(points); m_NumberOfConnectionsAdded = 0; m_InterpolationResult = nullptr; m_CurrentNumberOfReducedContours = 0; } mitk::SurfaceInterpolationController::~SurfaceInterpolationController() { // Removing all observers this->RemoveObservers(); } void mitk::SurfaceInterpolationController::RemoveObservers() { // Removing all observers auto dataIter = m_SegmentationObserverTags.begin(); for (; dataIter != m_SegmentationObserverTags.end(); ++dataIter) { (*dataIter).first->RemoveObserver((*dataIter).second); } m_SegmentationObserverTags.clear(); } mitk::SurfaceInterpolationController *mitk::SurfaceInterpolationController::GetInstance() { static mitk::SurfaceInterpolationController::Pointer m_Instance; if (m_Instance.IsNull()) { m_Instance = SurfaceInterpolationController::New(); } return m_Instance; } void mitk::SurfaceInterpolationController::AddNewContour(mitk::Surface::Pointer newContour) { if (newContour->GetVtkPolyData()->GetNumberOfPoints() > 0) { ContourPositionInformation contourInfo = CreateContourPositionInformation(newContour, nullptr); this->AddToInterpolationPipeline(contourInfo); this->Modified(); } } void mitk::SurfaceInterpolationController::AddNewContours(const std::vector& newContours, std::vector& contourPlanes, bool reinitializationAction) { if (nullptr == m_SelectedSegmentation) return; if (newContours.size() != contourPlanes.size()) { MITK_ERROR << "SurfaceInterpolationController::AddNewContours. contourPlanes and newContours are not of the same size."; } for (size_t i = 0; i < newContours.size(); ++i) { const auto &newContour = newContours[i]; const mitk::PlaneGeometry * planeGeometry = contourPlanes[i]; if (newContour->GetVtkPolyData()->GetNumberOfPoints() > 0) { auto contourInfo = CreateContourPositionInformation(newContour, planeGeometry); if (!reinitializationAction) { contourInfo.ContourPoint = this->ComputeInteriorPointOfContour(contourInfo, dynamic_cast(m_SelectedSegmentation) ); } else { auto vtkPolyData = contourInfo.Contour->GetVtkPolyData(); auto pointVtkArray = vtkDoubleArray::SafeDownCast(vtkPolyData->GetFieldData()->GetAbstractArray(1)); mitk::ScalarType *ptArr = new mitk::ScalarType[3]; for (int i = 0; i < pointVtkArray->GetSize(); ++i) ptArr[i] = pointVtkArray->GetValue(i); mitk::Point3D pt3D; pt3D.FillPoint(ptArr); contourInfo.ContourPoint = pt3D; } this->AddToInterpolationPipeline(contourInfo, reinitializationAction); } } this->Modified(); } mitk::DataNode* mitk::SurfaceInterpolationController::GetSegmentationImageNode() { DataNode* segmentationNode = nullptr; mitk::NodePredicateDataUID::Pointer dataUIDPredicate = mitk::NodePredicateDataUID::New(m_SelectedSegmentation->GetUID()); auto dataNodeObjects = m_DataStorage->GetSubset(dataUIDPredicate); if (dataNodeObjects->Size() != 0) { for (auto it = dataNodeObjects->Begin(); it != dataNodeObjects->End(); ++it) { segmentationNode = it->Value(); } } else { MITK_ERROR << "Unable to find the labelSetImage with the desired UID."; } return segmentationNode; } void mitk::SurfaceInterpolationController::AddPlaneGeometryNodeToDataStorage(const ContourPositionInformation& contourInfo) { auto planeGeometry = contourInfo.Plane; auto planeGeometryData = mitk::PlanarCircle::New(); planeGeometryData->SetPlaneGeometry(planeGeometry); mitk::Point2D p1; planeGeometry->Map(planeGeometry->GetCenter(), p1); planeGeometryData->PlaceFigure(p1); planeGeometryData->SetCurrentControlPoint(p1); planeGeometryData->SetProperty("initiallyplaced", mitk::BoolProperty::New(true)); if (planeGeometry) { auto segmentationNode = this->GetSegmentationImageNode(); auto isContourPlaneGeometry = mitk::NodePredicateProperty::New("isContourPlaneGeometry", mitk::BoolProperty::New(true)); mitk::DataStorage::SetOfObjects::ConstPointer contourNodes = m_DataStorage->GetDerivations(segmentationNode, isContourPlaneGeometry); auto contourFound = false; // Go through the pre-existing contours and check if the contour position matches them. for (auto it = contourNodes->Begin(); it != contourNodes->End(); ++it) { auto layerID = dynamic_cast(it->Value()->GetProperty("layerID"))->GetValue(); auto labelID = dynamic_cast(it->Value()->GetProperty("labelID"))->GetValue(); auto posID = dynamic_cast(it->Value()->GetProperty("position"))->GetValue(); bool sameLayer = (layerID == contourInfo.LayerValue); bool sameLabel = (labelID == contourInfo.LabelValue); bool samePos = (posID == contourInfo.Pos); if (samePos & sameLabel & sameLayer) { contourFound = true; it->Value()->SetData(planeGeometryData); break; } } if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { MITK_ERROR << "Invalid time point requested in AddPlaneGeometryNodeToDataStorage."; return; } const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); // Go through the contourPlaneGeometry Data and add the segmentationNode to it. if (!contourFound) { std::string contourName = "contourPlane " + std::to_string(m_ContourIndex); auto contourPlaneGeometryDataNode = mitk::DataNode::New(); contourPlaneGeometryDataNode->SetData(planeGeometryData); // No need to change properties contourPlaneGeometryDataNode->SetProperty("helper object", mitk::BoolProperty::New(false)); contourPlaneGeometryDataNode->SetProperty("hidden object", mitk::BoolProperty::New(true)); contourPlaneGeometryDataNode->SetProperty("isContourPlaneGeometry", mitk::BoolProperty::New(true)); contourPlaneGeometryDataNode->SetVisibility(false); // Need to change properties contourPlaneGeometryDataNode->SetProperty("name", mitk::StringProperty::New(contourName) ); contourPlaneGeometryDataNode->SetProperty("layerID", mitk::UIntProperty::New(contourInfo.LayerValue)); contourPlaneGeometryDataNode->SetProperty("labelID", mitk::UShortProperty::New(contourInfo.LabelValue)); contourPlaneGeometryDataNode->SetProperty("position", mitk::IntProperty::New(contourInfo.Pos)); contourPlaneGeometryDataNode->SetProperty("timeStep", mitk::IntProperty::New(currentTimeStep)); contourPlaneGeometryDataNode->SetProperty("px", mitk::DoubleProperty::New(contourInfo.ContourPoint[0])); contourPlaneGeometryDataNode->SetProperty("py", mitk::DoubleProperty::New(contourInfo.ContourPoint[1])); contourPlaneGeometryDataNode->SetProperty("pz", mitk::DoubleProperty::New(contourInfo.ContourPoint[2])); m_DataStorage->Add(contourPlaneGeometryDataNode, segmentationNode); } } } void mitk::SurfaceInterpolationController::AddToInterpolationPipeline(ContourPositionInformation& contourInfo, bool reinitializationAction) { if (!m_SelectedSegmentation) return; if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { MITK_ERROR << "Invalid time point requested for interpolation pipeline."; return; } // Get current time step either from the auto GetCurrentTimeStep = [=](ContourPositionInformation contourInfo) { if (reinitializationAction) { return contourInfo.TimeStep; } return static_cast(m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint)); }; const auto currentTimeStep = GetCurrentTimeStep(contourInfo); auto GetContourLayerID = [=](ContourPositionInformation contourInfo) { unsigned int currentLayerID; if(reinitializationAction) { if (contourInfo.LayerValue == std::numeric_limits::max()) { MITK_ERROR << "In mitk::SurfaceInterpolationController::AddToInterpolationPipeline. Problem in finding layerID"; } currentLayerID = contourInfo.LayerValue; } else { try { currentLayerID = dynamic_cast(m_SelectedSegmentation)->GetActiveLayer(); } catch (const std::exception& e) { MITK_ERROR << "Unable to cast image to LabelSetImage. " << e.what() << '\n'; } } return currentLayerID; }; unsigned int currentLayerID = GetContourLayerID(contourInfo); ContourPositionInformationVec3D ¤tImageContours = m_ListOfContours.at(m_SelectedSegmentation); ContourPositionInformationVec2D ¤tTimeStepContoursList = currentImageContours.at(currentTimeStep); ContourPositionInformationList ¤tContourList = currentTimeStepContoursList.at(currentLayerID); int replacementIndex = -1; int pos = -1; mitk::Surface* newContour = contourInfo.Contour; for (size_t i = 0; i < currentContourList.size(); i++) { auto& contourFromList = currentContourList.at(i); bool contoursAreCoplanar = ContoursCoplanar(contourInfo, contourFromList); bool contoursHaveSameLabel = contourInfo.LabelValue == contourFromList.LabelValue; // Coplanar contours have the same "pos". if (contoursAreCoplanar) { pos = contourFromList.Pos; if (contoursHaveSameLabel) { replacementIndex = i; } } } // The current contour has the same label and position as the current slice and a replacement is done. if (replacementIndex != -1) { contourInfo.Pos = pos; m_ListOfContours.at(m_SelectedSegmentation).at(currentTimeStep).at(currentLayerID).at(replacementIndex) = contourInfo; if (!reinitializationAction) { this->AddPlaneGeometryNodeToDataStorage(contourInfo); } return; } // Case that there is no contour in the current slice with the current label if (pos == -1) pos = m_ContourPosIndex++; m_ContourIndex++; contourInfo.Pos = pos; m_ListOfContours.at(m_SelectedSegmentation).at(currentTimeStep).at(currentLayerID).push_back(contourInfo); if (contourInfo.Plane == nullptr) { MITK_ERROR << "contourInfo plane is null."; } if (!reinitializationAction) { this->AddPlaneGeometryNodeToDataStorage(contourInfo); } if (newContour->GetVtkPolyData()->GetNumberOfPoints() == 0) { this->RemoveContour(contourInfo); if (m_ContourIndex > 0) m_ContourIndex--; if (m_ContourIndex > 0) m_ContourIndex--; } } bool mitk::SurfaceInterpolationController::RemoveContour(ContourPositionInformation contourInfo) { if (!m_SelectedSegmentation) { return false; } if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { return false; } const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); unsigned int currentLayerID = 0; try { currentLayerID = dynamic_cast(m_SelectedSegmentation)->GetActiveLayer(); } catch (const std::exception& e) { MITK_ERROR << e.what() << '\n'; } auto it = m_ListOfContours.at(m_SelectedSegmentation).at(currentTimeStep).at(currentLayerID).begin(); while (it != m_ListOfContours.at(m_SelectedSegmentation).at(currentTimeStep).at(currentLayerID).end()) { const ContourPositionInformation ¤tContour = (*it); if (ContoursCoplanar(currentContour, contourInfo)) { m_ListOfContours.at(m_SelectedSegmentation).at(currentTimeStep).at(currentLayerID).erase(it); this->ReinitializeInterpolation(); return true; } ++it; } return false; } const mitk::Surface *mitk::SurfaceInterpolationController::GetContour(const ContourPositionInformation &contourInfo) { if (!m_SelectedSegmentation) { return nullptr; } if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { return nullptr; } const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); const auto activeLayerID = dynamic_cast(m_SelectedSegmentation)->GetActiveLayer(); const auto &contourList = m_ListOfContours.at(m_SelectedSegmentation).at(currentTimeStep).at(activeLayerID); for (auto ¤tContour : contourList) { if (ContoursCoplanar(contourInfo, currentContour)) { return currentContour.Contour; } } return nullptr; } unsigned int mitk::SurfaceInterpolationController::GetNumberOfContours() { if (!m_SelectedSegmentation) { return -1; } if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { return -1; } const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); auto contourDoubleList = m_ListOfContours.at(m_SelectedSegmentation).at(currentTimeStep); unsigned int numContours = 0; for (auto& contourList : contourDoubleList) { numContours += contourList.size(); } return numContours; } void mitk::SurfaceInterpolationController::AddActiveLabelContoursForInterpolation(mitk::Label::PixelType activeLabel) { this->ReinitializeInterpolation(); if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { MITK_ERROR << "Invalid time point requested for interpolation pipeline."; return; } const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); unsigned int currentLayerID = 0; try { currentLayerID = dynamic_cast(m_SelectedSegmentation)->GetActiveLayer(); } catch (const std::exception& e) { MITK_ERROR << e.what() << '\n'; } ContourPositionInformationVec3D ¤tImageContours = m_ListOfContours.at(m_SelectedSegmentation); if (currentImageContours.size() <= currentTimeStep) { MITK_INFO << "Contours for current time step don't exist."; return; } ContourPositionInformationVec2D ¤tTimeStepContoursList = currentImageContours.at(currentTimeStep); if (currentTimeStepContoursList.size() <= currentLayerID) { MITK_INFO << "Contours for current layer don't exist."; return; } ContourPositionInformationList ¤tContours = currentTimeStepContoursList.at(currentLayerID); for (size_t i = 0; i < currentContours.size(); ++i) { if (currentContours.at(i).LabelValue == activeLabel) { m_ListOfInterpolationSessions.at(m_SelectedSegmentation).at(currentTimeStep).push_back(currentContours.at(i)); m_ReduceFilter->SetInput(m_ListOfInterpolationSessions.at(m_SelectedSegmentation).at(currentTimeStep).size()-1, currentContours.at(i).Contour); } } } void mitk::SurfaceInterpolationController::Interpolate() { if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { MITK_WARN << "No interpolation possible, currently selected timepoint is not in the time bounds of currently selected segmentation. Time point: " << m_CurrentTimePoint; m_InterpolationResult = nullptr; return; } const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); m_ReduceFilter->Update(); m_CurrentNumberOfReducedContours = m_ReduceFilter->GetNumberOfOutputs(); if (m_CurrentNumberOfReducedContours == 1) { vtkPolyData *tmp = m_ReduceFilter->GetOutput(0)->GetVtkPolyData(); if (tmp == nullptr) { m_CurrentNumberOfReducedContours = 0; } } // We use the timeSelector to get the segmentation image for the current segmentation. mitk::ImageTimeSelector::Pointer timeSelector = mitk::ImageTimeSelector::New(); timeSelector->SetInput(m_SelectedSegmentation); timeSelector->SetTimeNr(currentTimeStep); timeSelector->SetChannelNr(0); timeSelector->Update(); mitk::Image::Pointer refSegImage = timeSelector->GetOutput(); itk::ImageBase<3>::Pointer itkImage = itk::ImageBase<3>::New(); AccessFixedDimensionByItk_1(refSegImage, GetImageBase, 3, itkImage); m_NormalsFilter->SetSegmentationBinaryImage(refSegImage); for (size_t i = 0; i < m_CurrentNumberOfReducedContours; ++i) { mitk::Surface::Pointer reducedContour = m_ReduceFilter->GetOutput(i); reducedContour->DisconnectPipeline(); m_NormalsFilter->SetInput(i, reducedContour); m_InterpolateSurfaceFilter->SetInput(i, m_NormalsFilter->GetOutput(i)); } if (m_CurrentNumberOfReducedContours < 2) { // If no interpolation is possible reset the interpolation result MITK_INFO << "Interpolation impossible: not enough contours."; m_InterpolationResult = nullptr; return; } // Setting up progress bar mitk::ProgressBar::GetInstance()->AddStepsToDo(10); // create a surface from the distance-image mitk::ImageToSurfaceFilter::Pointer imageToSurfaceFilter = mitk::ImageToSurfaceFilter::New(); imageToSurfaceFilter->SetInput(m_InterpolateSurfaceFilter->GetOutput()); imageToSurfaceFilter->SetThreshold(0); imageToSurfaceFilter->SetSmooth(true); imageToSurfaceFilter->SetSmoothIteration(1); imageToSurfaceFilter->Update(); mitk::Surface::Pointer interpolationResult = mitk::Surface::New(); interpolationResult->Expand(m_SelectedSegmentation->GetTimeSteps()); auto geometry = m_SelectedSegmentation->GetTimeGeometry()->Clone(); geometry->ReplaceTimeStepGeometries(mitk::Geometry3D::New()); interpolationResult->SetTimeGeometry(geometry); interpolationResult->SetVtkPolyData(imageToSurfaceFilter->GetOutput()->GetVtkPolyData(), currentTimeStep); m_InterpolationResult = interpolationResult; m_DistanceImageSpacing = m_InterpolateSurfaceFilter->GetDistanceImageSpacing(); auto* contoursGeometry = static_cast(m_Contours->GetTimeGeometry()); auto timeBounds = geometry->GetTimeBounds(currentTimeStep); contoursGeometry->SetFirstTimePoint(timeBounds[0]); contoursGeometry->SetStepDuration(timeBounds[1] - timeBounds[0]); // Last progress step mitk::ProgressBar::GetInstance()->Progress(20); m_InterpolationResult->DisconnectPipeline(); } mitk::Surface::Pointer mitk::SurfaceInterpolationController::GetInterpolationResult() { return m_InterpolationResult; } mitk::Surface *mitk::SurfaceInterpolationController::GetContoursAsSurface() { return m_Contours; } void mitk::SurfaceInterpolationController::SetDataStorage(DataStorage::Pointer ds) { m_DataStorage = ds; } void mitk::SurfaceInterpolationController::SetMinSpacing(double minSpacing) { m_ReduceFilter->SetMinSpacing(minSpacing); } void mitk::SurfaceInterpolationController::SetMaxSpacing(double maxSpacing) { m_ReduceFilter->SetMaxSpacing(maxSpacing); m_NormalsFilter->SetMaxSpacing(maxSpacing); } void mitk::SurfaceInterpolationController::SetDistanceImageVolume(unsigned int distImgVolume) { m_InterpolateSurfaceFilter->SetDistanceImageVolume(distImgVolume); } mitk::Image::Pointer mitk::SurfaceInterpolationController::GetCurrentSegmentation() { return m_SelectedSegmentation; } mitk::Image *mitk::SurfaceInterpolationController::GetImage() { return m_InterpolateSurfaceFilter->GetOutput(); } double mitk::SurfaceInterpolationController::EstimatePortionOfNeededMemory() { double numberOfPointsAfterReduction = m_ReduceFilter->GetNumberOfPointsAfterReduction() * 3; double sizeOfPoints = pow(numberOfPointsAfterReduction, 2) * sizeof(double); double totalMem = mitk::MemoryUtilities::GetTotalSizeOfPhysicalRam(); double percentage = sizeOfPoints / totalMem; return percentage; } unsigned int mitk::SurfaceInterpolationController::GetNumberOfInterpolationSessions() { return m_ListOfInterpolationSessions.size(); } template void mitk::SurfaceInterpolationController::GetImageBase(itk::Image *input, itk::ImageBase<3>::Pointer &result) { result->Graft(input); } void mitk::SurfaceInterpolationController::SetCurrentSegmentationInterpolationList(mitk::Image::Pointer segmentation) { this->SetCurrentInterpolationSession(segmentation); } void mitk::SurfaceInterpolationController::SetCurrentInterpolationSession(mitk::Image::Pointer currentSegmentationImage) { if (currentSegmentationImage.GetPointer() == m_SelectedSegmentation) { return; } if (currentSegmentationImage.IsNull()) { m_SelectedSegmentation = nullptr; return; } m_SelectedSegmentation = currentSegmentationImage.GetPointer(); try { auto labelSetImage = dynamic_cast(m_SelectedSegmentation); auto it = m_ListOfContours.find(currentSegmentationImage.GetPointer()); // If the session does not exist yet create a new ContourPositionPairList otherwise reinitialize the interpolation // pipeline if (it == m_ListOfContours.end()) { ContourPositionInformationVec3D newList; auto numTimeSteps = labelSetImage->GetTimeGeometry()->CountTimeSteps(); for (size_t t = 0; t < numTimeSteps; ++t) { auto twoDList = ContourPositionInformationVec2D(); auto contourList = ContourPositionInformationList(); twoDList.push_back(contourList); newList.push_back(twoDList); } m_ListOfContours[m_SelectedSegmentation] = newList; m_InterpolationResult = nullptr; m_CurrentNumberOfReducedContours = 0; auto command = itk::MemberCommand::New(); command->SetCallbackFunction(this, &SurfaceInterpolationController::OnSegmentationDeleted); m_SegmentationObserverTags[m_SelectedSegmentation] = labelSetImage->AddObserver(itk::DeleteEvent(), command); m_NumberOfLayersInCurrentSegmentation = labelSetImage->GetNumberOfLayers(); } // auto labelSetImage = dynamic_cast(m_SelectedSegmentation); auto numLayersInSelectedSegmentation = labelSetImage->GetNumberOfLayers(); // Maybe this has to change. for (size_t layerID = 0; layerID < numLayersInSelectedSegmentation; ++layerID) { this->AddLabelSetConnection(layerID); } } catch (const std::exception &e) { MITK_ERROR << "Unable to cast image as LabelSetImage"; } auto it2 = m_ListOfInterpolationSessions.find(currentSegmentationImage.GetPointer()); if (it2 == m_ListOfInterpolationSessions.end()) { ContourPositionInformationVec2D newList; m_ListOfInterpolationSessions[m_SelectedSegmentation] = newList; m_InterpolationResult = nullptr; m_CurrentNumberOfReducedContours = 0; } this->ReinitializeInterpolation(); } bool mitk::SurfaceInterpolationController::ReplaceInterpolationSession(mitk::Image::Pointer oldSession, mitk::Image::Pointer newSession) { if (oldSession.IsNull() || newSession.IsNull()) return false; if (oldSession.GetPointer() == newSession.GetPointer()) return false; if (!mitk::Equal(*(oldSession->GetGeometry()), *(newSession->GetGeometry()), mitk::eps, false)) return false; auto it = m_ListOfInterpolationSessions.find(oldSession.GetPointer()); if (it == m_ListOfInterpolationSessions.end()) return false; if (!newSession->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { MITK_WARN << "Interpolation session cannot be replaced. Currently selected timepoint is not in the time bounds of the new session. Time point: " << m_CurrentTimePoint; return false; } ContourPositionInformationVec2D oldList = (*it).second; m_ListOfInterpolationSessions[newSession.GetPointer()] = oldList; itk::MemberCommand::Pointer command = itk::MemberCommand::New(); command->SetCallbackFunction(this, &SurfaceInterpolationController::OnSegmentationDeleted); m_SegmentationObserverTags[newSession] = newSession->AddObserver(itk::DeleteEvent(), command); if (m_SelectedSegmentation == oldSession) m_SelectedSegmentation = newSession; const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); mitk::ImageTimeSelector::Pointer timeSelector = mitk::ImageTimeSelector::New(); timeSelector->SetInput(m_SelectedSegmentation); timeSelector->SetTimeNr(currentTimeStep); timeSelector->SetChannelNr(0); timeSelector->Update(); mitk::Image::Pointer refSegImage = timeSelector->GetOutput(); m_NormalsFilter->SetSegmentationBinaryImage(refSegImage); this->RemoveInterpolationSession(oldSession); return true; } void mitk::SurfaceInterpolationController::RemoveSegmentationFromContourList(mitk::Image *segmentation) { this->RemoveInterpolationSession(segmentation); } void mitk::SurfaceInterpolationController::RemoveInterpolationSession(mitk::Image::Pointer segmentationImage) { if (segmentationImage) { if (m_SelectedSegmentation == segmentationImage) { m_NormalsFilter->SetSegmentationBinaryImage(nullptr); m_SelectedSegmentation = nullptr; } m_ListOfInterpolationSessions.erase(segmentationImage); m_ListOfContours.erase(segmentationImage); // Remove observer auto pos = m_SegmentationObserverTags.find(segmentationImage); if (pos != m_SegmentationObserverTags.end()) { segmentationImage->RemoveObserver((*pos).second); m_SegmentationObserverTags.erase(pos); } } } void mitk::SurfaceInterpolationController::RemoveAllInterpolationSessions() { // Removing all observers auto dataIter = m_SegmentationObserverTags.begin(); while (dataIter != m_SegmentationObserverTags.end()) { mitk::Image *image = (*dataIter).first; image->RemoveObserver((*dataIter).second); ++dataIter; } m_SegmentationObserverTags.clear(); m_SelectedSegmentation = nullptr; m_ListOfInterpolationSessions.clear(); m_ListOfContours.clear(); } template std::vector GetPixelValuesPresentInImage(mitk::LabelSetImage* labelSetImage) { mitk::ImagePixelReadAccessor readAccessor(labelSetImage); std::vector pixelsPresent; std::size_t numberOfPixels = 1; for (int dim = 0; dim < static_cast(VImageDimension); ++dim) numberOfPixels *= static_cast(readAccessor.GetDimension(dim)); auto src = readAccessor.GetData(); for (std::size_t i = 0; i < numberOfPixels; ++i) { mitk::Label::PixelType pixelVal = *(src + i); if ( (std::find(pixelsPresent.begin(), pixelsPresent.end(), pixelVal) == pixelsPresent.end()) && (pixelVal != 0) ) { pixelsPresent.push_back(pixelVal); } } return pixelsPresent; } void mitk::SurfaceInterpolationController::RemoveContours(mitk::Label::PixelType label, unsigned int timeStep, unsigned int layerID) { auto isContourEqualToLabelValue = [label] (ContourPositionInformation& contour) -> bool { return (contour.LabelValue == label); }; ContourPositionInformationVec3D ¤tImageContours = m_ListOfContours.at(m_SelectedSegmentation); ContourPositionInformationList ¤tContourList = currentImageContours.at(timeStep).at(layerID); unsigned int numContoursBefore = currentContourList.size(); auto it = std::remove_if(currentContourList.begin(), currentContourList.end(), isContourEqualToLabelValue); currentContourList.erase(it, currentContourList.end()); unsigned int numContoursAfter = currentContourList.size(); unsigned int numContours = numContoursAfter - numContoursBefore; m_ContourIndex -= numContours; } void mitk::SurfaceInterpolationController::OnSegmentationDeleted(const itk::Object *caller, const itk::EventObject & /*event*/) { auto *tempImage = dynamic_cast(const_cast(caller)); if (tempImage) { if (m_SelectedSegmentation == tempImage) { m_NormalsFilter->SetSegmentationBinaryImage(nullptr); m_SelectedSegmentation = nullptr; } m_SegmentationObserverTags.erase(tempImage); m_ListOfContours.erase(tempImage); m_ListOfInterpolationSessions.erase(tempImage); } } void mitk::SurfaceInterpolationController::ReinitializeInterpolation() { // If session has changed reset the pipeline m_ReduceFilter->Reset(); m_NormalsFilter->Reset(); m_InterpolateSurfaceFilter->Reset(); // Empty out the listOfInterpolationSessions m_ListOfInterpolationSessions[m_SelectedSegmentation].clear(); itk::ImageBase<3>::Pointer itkImage = itk::ImageBase<3>::New(); if (m_SelectedSegmentation) { if (!m_SelectedSegmentation->GetTimeGeometry()->IsValidTimePoint(m_CurrentTimePoint)) { MITK_WARN << "Interpolation cannot be reinitialized. Currently selected timepoint is not in the time bounds of the currently selected segmentation. Time point: " << m_CurrentTimePoint; return; } const auto currentTimeStep = m_SelectedSegmentation->GetTimeGeometry()->TimePointToTimeStep(m_CurrentTimePoint); // Set reference image for interpolation surface filter mitk::ImageTimeSelector::Pointer timeSelector = mitk::ImageTimeSelector::New(); timeSelector->SetInput(m_SelectedSegmentation); timeSelector->SetTimeNr(currentTimeStep); timeSelector->SetChannelNr(0); timeSelector->Update(); mitk::Image::Pointer refSegImage = timeSelector->GetOutput(); AccessFixedDimensionByItk_1(refSegImage, GetImageBase, 3, itkImage); m_InterpolateSurfaceFilter->SetReferenceImage(itkImage.GetPointer()); // Resize listofinterpolationsessions and listofcontours to numTimeSteps unsigned int numTimeSteps = m_SelectedSegmentation->GetTimeSteps(); unsigned int size = m_ListOfInterpolationSessions[m_SelectedSegmentation].size(); if (size != numTimeSteps) { m_ListOfInterpolationSessions.at(m_SelectedSegmentation).resize(numTimeSteps); } } } void mitk::SurfaceInterpolationController::AddLabelSetConnection(unsigned int layerID) { if (m_SelectedSegmentation != nullptr) { try { auto workingImage = dynamic_cast(m_SelectedSegmentation); auto previousLayerID = workingImage->GetActiveLayer(); workingImage->SetActiveLayer(layerID); auto activeLabelSet = workingImage->GetLabelSet(layerID); activeLabelSet->RemoveLabelEvent += mitk::MessageDelegate1( this, &mitk::SurfaceInterpolationController::OnRemoveLabel); activeLabelSet->ActiveLabelEvent += mitk::MessageDelegate1( this, &mitk::SurfaceInterpolationController::OnActiveLabel); workingImage->AfterChangeLayerEvent += mitk::MessageDelegate( this, &mitk::SurfaceInterpolationController::OnLayerChanged); m_NumberOfConnectionsAdded += 1; workingImage->SetActiveLayer(previousLayerID); } catch(const std::exception& e) { MITK_ERROR << e.what() << '\n'; } } } void mitk::SurfaceInterpolationController::AddLabelSetConnection() { if (m_SelectedSegmentation != nullptr) { try { auto workingImage = dynamic_cast(m_SelectedSegmentation); auto activeLabelSet = workingImage->GetActiveLabelSet(); activeLabelSet->RemoveLabelEvent += mitk::MessageDelegate1( this, &mitk::SurfaceInterpolationController::OnRemoveLabel); workingImage->GetActiveLabelSet()->ActiveLabelEvent += mitk::MessageDelegate1( this, &mitk::SurfaceInterpolationController::OnActiveLabel); workingImage->AfterChangeLayerEvent += mitk::MessageDelegate( this, &mitk::SurfaceInterpolationController::OnLayerChanged); m_NumberOfConnectionsAdded += 1; } catch(const std::exception& e) { MITK_ERROR << e.what() << '\n'; } } } void mitk::SurfaceInterpolationController::RemoveLabelSetConnection(mitk::LabelSetImage* labelSetImage, unsigned int layerID) { labelSetImage->SetActiveLayer(layerID); labelSetImage->GetActiveLabelSet()->RemoveLabelEvent -= mitk::MessageDelegate1( this, &mitk::SurfaceInterpolationController::OnRemoveLabel); // labelSetImage->GetActiveLabelSet()->ActiveLabelEvent -= mitk::MessageDelegate1( // this, &mitk::SurfaceInterpolationController::OnActiveLabel); labelSetImage->AfterChangeLayerEvent -= mitk::MessageDelegate( this, &mitk::SurfaceInterpolationController::OnLayerChanged); m_NumberOfConnectionsAdded -= 1; } void mitk::SurfaceInterpolationController::RemoveLabelSetConnection() { if (m_SelectedSegmentation != nullptr) { try { auto workingImage = dynamic_cast(m_SelectedSegmentation); workingImage->GetActiveLabelSet()->RemoveLabelEvent -= mitk::MessageDelegate1( this, &mitk::SurfaceInterpolationController::OnRemoveLabel); workingImage->GetActiveLabelSet()->ActiveLabelEvent -= mitk::MessageDelegate1( this, &mitk::SurfaceInterpolationController::OnActiveLabel); workingImage->AfterChangeLayerEvent -= mitk::MessageDelegate( this, &mitk::SurfaceInterpolationController::OnLayerChanged); } catch (const std::exception& e) { std::cerr << e.what() << '\n'; } } } void mitk::SurfaceInterpolationController::OnRemoveLabel(mitk::Label::PixelType /*removedLabelValue*/) { if (m_SelectedSegmentation != nullptr) { auto numTimeSteps = m_SelectedSegmentation->GetTimeGeometry()->CountTimeSteps(); try { auto labelSetImage = dynamic_cast(m_SelectedSegmentation); auto currentLayerID = labelSetImage->GetActiveLayer(); for(unsigned int t = 0; t < numTimeSteps; ++t) { this->RemoveContours(m_PreviousActiveLabelValue,t,currentLayerID); } } catch(const std::exception& e) { std::cerr << e.what() << '\n'; } } } void mitk::SurfaceInterpolationController::OnActiveLabel(mitk::Label::PixelType newActiveLabelValue) { m_PreviousActiveLabelValue = m_CurrentActiveLabelValue; m_CurrentActiveLabelValue = newActiveLabelValue; } unsigned int mitk::SurfaceInterpolationController::GetNumberOfLayersInCurrentSegmentation() const { return m_NumberOfLayersInCurrentSegmentation; } void mitk::SurfaceInterpolationController::SetNumberOfLayersInCurrentSegmentation(unsigned int numLayers) { m_NumberOfLayersInCurrentSegmentation = numLayers; } void mitk::SurfaceInterpolationController::OnAddLayer() { assert(m_SelectedSegmentation != nullptr); auto& contoursForSegmentation = m_ListOfContours.at(m_SelectedSegmentation); // Push an information list for each time step. for(size_t t = 0; t < contoursForSegmentation.size(); ++t) { contoursForSegmentation.at(t).push_back( ContourPositionInformationList() ); } } void mitk::SurfaceInterpolationController::OnRemoveLayer() { assert(m_SelectedSegmentation != nullptr); auto& contoursForSegmentation = m_ListOfContours.at(m_SelectedSegmentation); // Erase the layers in each of the time steps. // The previous layer is removed for (size_t t = 0; t < contoursForSegmentation.size(); ++t) { assert(m_PreviousLayerIndex < contoursForSegmentation.at(t).size()); auto& contoursAtTimeStep = contoursForSegmentation.at(t); for (size_t c = m_CurrentLayerIndex+1; c < contoursAtTimeStep.size(); ++c) { auto& contoursInCurrentLayer = contoursAtTimeStep.at(c); for (auto& contour : contoursInCurrentLayer) { contour.LayerValue = contour.LayerValue - 1; } } } for (size_t t = 0; t < contoursForSegmentation.size(); ++t) { assert (m_CurrentLayerIndex < contoursForSegmentation.at(t).size()); contoursForSegmentation.at(t).erase(contoursForSegmentation.at(t).begin() + m_PreviousLayerIndex); } this->Modified(); } void mitk::SurfaceInterpolationController::OnLayerChanged() { auto currentLayer = dynamic_cast(m_SelectedSegmentation)->GetActiveLayer(); m_PreviousLayerIndex = m_CurrentLayerIndex; m_CurrentLayerIndex = currentLayer; } -mitk::SurfaceInterpolationController::ContourPositionInformationList& mitk::SurfaceInterpolationController::GetContours(unsigned int timeStep, unsigned int layerID) +mitk::SurfaceInterpolationController::ContourPositionInformationList* mitk::SurfaceInterpolationController::GetContours(unsigned int timeStep, unsigned int layerID) { if (m_SelectedSegmentation == nullptr) - { - MITK_ERROR << "Invalid segmentation from mitk::SurfaceInterpolationController::GetContours"; - } + return nullptr; + if (timeStep >= m_ListOfContours.at(m_SelectedSegmentation).size()) - { - MITK_ERROR << "Invalid timeStep from mitk::SurfaceInterpolationController::GetContours"; - } + return nullptr; + if (layerID >= m_ListOfContours.at(m_SelectedSegmentation).at(timeStep).size()) - { - MITK_ERROR << "Invalid timeStep from mitk::SurfaceInterpolationController::GetContours"; - } - return m_ListOfContours.at(m_SelectedSegmentation).at(timeStep).at(layerID); + return nullptr; + + return &m_ListOfContours[m_SelectedSegmentation][timeStep][layerID]; } void mitk::SurfaceInterpolationController::CompleteReinitialization(const std::vector& contourList, std::vector& contourPlanes) { this->ClearInterpolationSession(); auto labelSetImage = dynamic_cast(m_SelectedSegmentation); auto numLayers = labelSetImage->GetNumberOfLayers(); // Add layers to the m_ListOfContours for (size_t layer = 0; layer < numLayers; ++layer) { this->OnAddLayer(); } // Now the layers should be empty and the new layers can be added. this->AddNewContours(contourList, contourPlanes, true); } void mitk::SurfaceInterpolationController::ClearInterpolationSession() { if (m_SelectedSegmentation != nullptr) { auto it = m_ListOfContours.find(m_SelectedSegmentation); if (it != m_ListOfContours.end()) { auto timeSteps = m_ListOfContours[m_SelectedSegmentation].size(); try { auto labelSetImage = dynamic_cast(m_SelectedSegmentation); auto labelSetImageTimeSteps = labelSetImage->GetTimeGeometry()->CountTimeSteps(); if (timeSteps != labelSetImageTimeSteps) { MITK_ERROR << "Time steps are not the same."; } for (size_t t = 0; t < timeSteps; ++t) { m_ListOfContours[m_SelectedSegmentation][t].clear(); } } catch(std::bad_cast& e) { MITK_ERROR << "Unable to cast m_SelectedSegmentation to labelSetImage in ClearInterpolationSession"; } } } } std::vector< mitk::Point3D > mitk::ContourExt::GetBoundingBoxGridPoints( size_t planeDimension, double startDim1, size_t numPointsToSampleDim1, double deltaDim1, double startDim2, size_t numPointsToSampleDim2, double deltaDim2, double valuePlaneDim) { std::vector< mitk::Point3D > gridPoints; for (size_t i = 0; i < numPointsToSampleDim1; ++i) { for (size_t j = 0; j < numPointsToSampleDim2; ++j) { mitk::ScalarType *ptVec = new mitk::ScalarType[3]; if (planeDimension == 0) { ptVec[0] = valuePlaneDim; ptVec[1] = startDim1 + deltaDim1 * i; ptVec[2] = startDim2 + deltaDim2 * j; } else if (planeDimension == 1) { ptVec[0] = startDim1 + deltaDim1 * i; ptVec[1] = valuePlaneDim; ptVec[2] = startDim2 + deltaDim2 * j; } else if (planeDimension == 2) { ptVec[0] = startDim1 + deltaDim1 * i; ptVec[1] = startDim2 + deltaDim2 * j; ptVec[2] = valuePlaneDim; } mitk::Point3D pt3D; pt3D.FillPoint(ptVec); gridPoints.push_back(pt3D); } } return gridPoints; } mitk::Point3D mitk::SurfaceInterpolationController::ComputeInteriorPointOfContour( const mitk::SurfaceInterpolationController::ContourPositionInformation& contour, mitk::LabelSetImage * labelSetImage) { if (labelSetImage->GetDimension() == 4) { return mitk::ContourExt::ComputeInteriorPointOfContour<4>(contour, labelSetImage, m_CurrentTimePoint); } else { return mitk::ContourExt::ComputeInteriorPointOfContour<3>(contour, labelSetImage, m_CurrentTimePoint); } } template mitk::Point3D mitk::ContourExt::ComputeInteriorPointOfContour( const mitk::SurfaceInterpolationController::ContourPositionInformation& contour, mitk::LabelSetImage * labelSetImage, mitk::TimePointType currentTimePoint) { mitk::ImagePixelReadAccessor readAccessor(labelSetImage); if (!labelSetImage->GetTimeGeometry()->IsValidTimePoint(currentTimePoint)) { MITK_ERROR << "Invalid time point requested for interpolation pipeline."; mitk::Point3D pt; return pt; } std::vector pixelsPresent; const auto currentTimeStep = labelSetImage->GetTimeGeometry()->TimePointToTimeStep(currentTimePoint); auto polyData = contour.Contour->GetVtkPolyData(); polyData->ComputeCellsBounds(); mitk::ScalarType cellBounds[6]; polyData->GetCellsBounds(cellBounds); size_t numPointsToSample = 10; mitk::ScalarType StartX = cellBounds[0]; mitk::ScalarType StartY = cellBounds[2]; mitk::ScalarType StartZ = cellBounds[4]; size_t deltaX = (cellBounds[1] - cellBounds[0]) / numPointsToSample; size_t deltaY = (cellBounds[3] - cellBounds[2]) / numPointsToSample; size_t deltaZ = (cellBounds[5] - cellBounds[4]) / numPointsToSample; auto planeOrientation = mitk::ContourExt::GetContourOrientation(contour.ContourNormal); std::vector points; if (planeOrientation == 0) { points = mitk::ContourExt::GetBoundingBoxGridPoints(planeOrientation, StartY, numPointsToSample, deltaY, StartZ, numPointsToSample, deltaZ, StartX); } else if (planeOrientation == 1) { points = mitk::ContourExt::GetBoundingBoxGridPoints(planeOrientation, StartX, numPointsToSample, deltaX, StartZ, numPointsToSample, deltaZ, StartY); } else if (planeOrientation == 2) { points = mitk::ContourExt::GetBoundingBoxGridPoints(planeOrientation, StartX, numPointsToSample, deltaX, StartY, numPointsToSample, deltaY, StartZ); } mitk::Label::PixelType pixelVal; mitk::Point3D pt3D; std::vector pixelVals; for (size_t i = 0; i < points.size(); ++i) { pt3D = points[i]; itk::Index<3> itkIndex; labelSetImage->GetGeometry()->WorldToIndex(pt3D, itkIndex); if (VImageDimension == 4) { itk::Index time3DIndex; for (size_t i = 0; i < itkIndex.size(); ++i) time3DIndex[i] = itkIndex[i]; time3DIndex[3] = currentTimeStep; pixelVal = readAccessor.GetPixelByIndexSafe(time3DIndex); } else if (VImageDimension == 3) { itk::Index geomIndex; for (size_t i=0;i mitk::eps) { planeOrientation = 2; } else if (fabs(dotY) > mitk::eps) { planeOrientation = 1; } else if(fabs(dotX) > mitk::eps) { planeOrientation = 0; } return planeOrientation; } diff --git a/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.h b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.h index 5137c3a5ee..3a2ea6eccc 100644 --- a/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.h +++ b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.h @@ -1,468 +1,468 @@ /*============================================================================ The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center (DKFZ) All rights reserved. Use of this source code is governed by a 3-clause BSD license that can be found in the LICENSE file. ============================================================================*/ #ifndef mitkSurfaceInterpolationController_h #define mitkSurfaceInterpolationController_h #include #include #include #include #include namespace mitk { class ComputeContourSetNormalsFilter; class CreateDistanceImageFromSurfaceFilter; class LabelSetImage; class ReduceContourSetFilter; class MITKSURFACEINTERPOLATION_EXPORT SurfaceInterpolationController : public itk::Object { public: mitkClassMacroItkParent(SurfaceInterpolationController, itk::Object); itkFactorylessNewMacro(Self); itkCloneMacro(Self); itkGetMacro(DistanceImageSpacing, double); struct MITKSURFACEINTERPOLATION_EXPORT ContourPositionInformation { int Pos; unsigned int SliceIndex; Surface::Pointer Contour; Vector3D ContourNormal; Point3D ContourPoint; mitk::PlaneGeometry* Plane; mitk::Label::PixelType LabelValue; unsigned int LayerValue; size_t TimeStep; ContourPositionInformation() : Pos(-1), SliceIndex(0), Plane(nullptr), LabelValue(std::numeric_limits::max()), LayerValue(std::numeric_limits::max()), TimeStep(std::numeric_limits::max()) { } }; typedef std::vector ContourPositionInformationList; typedef std::vector ContourPositionInformationVec2D; // first index is the current time step. second index is the layerID. third index is the contour index. typedef std::vector ContourPositionInformationVec3D; typedef std::map ContourListMap; typedef std::map ContourContainer; static SurfaceInterpolationController *GetInstance(); void SetCurrentTimePoint(TimePointType tp) { if (m_CurrentTimePoint != tp) { m_CurrentTimePoint = tp; if (m_SelectedSegmentation) { this->ReinitializeInterpolation(); } } }; TimePointType GetCurrentTimePoint() const { return m_CurrentTimePoint; }; /** * @brief Adds a new extracted contour to the list * @param newContour the contour to be added. If a contour at that position * already exists the related contour will be updated */ void AddNewContour(Surface::Pointer newContour); /** * @brief Adds new extracted contours to the list. If one or more contours at a given position * already exist they will be updated respectively */ void AddNewContours(const std::vector& newContours, std::vector& contourPlanes, bool reinitializeAction = false); /** * @brief Returns the contour for a given plane for the current selected segmenation * @param contourInfo the contour which should be returned * @return the contour as an mitk::Surface. If no contour is available at the give position nullptr is returned */ const mitk::Surface *GetContour(const ContourPositionInformation& contourInfo); /** * @brief Computes an interior point of the input contour. It's used to detect merge and erase operations. * * @param contour Contour for which to compute the contour * @param labelSetImage LabelSetImage used input to check contour Label. * @return mitk::Point3D 3D Interior point of the contour returned. */ mitk::Point3D ComputeInteriorPointOfContour(const ContourPositionInformation& contour, mitk::LabelSetImage * labelSetImage); /** * @brief Make the surface interpolator responsive to the segmentation image by subscribing to events from the image. * */ void AddLabelSetConnection(); /** * @brief Make the surface interpolator responsive to the segmentation image by stopping subscription to events from the image. * */ void RemoveLabelSetConnection(); void RemoveLabelSetConnection(mitk::LabelSetImage* labelSetImage, unsigned int layerID); /** * @brief Resets the pipeline for interpolation. The various filters used are reset. * */ void ReinitializeInterpolation(); void RemoveObservers(); void AddLabelSetConnection(unsigned int layerID); void UnsetSelectedImage() { m_SelectedSegmentation = nullptr; } /** * @brief Returns the number of layers in the current segmentation image. * */ unsigned int GetNumberOfLayersInCurrentSegmentation() const; /** * @brief Set the number of layers in the current segmentation image. * */ void SetNumberOfLayersInCurrentSegmentation(unsigned int); /** * @brief Function that does the data management when a layer is removed. * */ void OnRemoveLayer(); /** * @brief Function that does the data management when a layer is added. * */ void OnAddLayer(); /** * @brief Returns the number of available contours for the current selected segmentation * @return the number of contours */ unsigned int GetNumberOfContours(); /** * @brief Performs the interpolation. * */ void Interpolate(); /** * @brief Get the Result of the interpolation operation. * * @return mitk::Surface::Pointer */ mitk::Surface::Pointer GetInterpolationResult(); /** * @brief Sets the minimum spacing of the current selected segmentation * This is needed since the contour points we reduced before they are used to interpolate the surface. * * @param minSpacing Paramter to set */ void SetMinSpacing(double minSpacing); /** * @brief Sets the minimum spacing of the current selected segmentation * This is needed since the contour points we reduced before they are used to interpolate the surface * @param maxSpacing Set the max Spacing for interpolation */ void SetMaxSpacing(double maxSpacing); /** * Sets the volume i.e. the number of pixels that the distance image should have * By evaluation we found out that 50.000 pixel delivers a good result */ void SetDistanceImageVolume(unsigned int distImageVolume); /** * @brief Get the current selected segmentation for which the interpolation is performed * @return the current segmentation image */ mitk::Image::Pointer GetCurrentSegmentation(); Surface *GetContoursAsSurface(); void SetDataStorage(DataStorage::Pointer ds); /** * Sets the current list of contourpoints which is used for the surface interpolation * @param segmentation The current selected segmentation * \deprecatedSince{2014_03} */ DEPRECATED(void SetCurrentSegmentationInterpolationList(mitk::Image::Pointer segmentation)); /** * Sets the current list of contourpoints which is used for the surface interpolation * @param currentSegmentationImage The current selected segmentation */ void SetCurrentInterpolationSession(mitk::Image::Pointer currentSegmentationImage); /** * Removes the segmentation and all its contours from the list * @param segmentation The segmentation to be removed * \deprecatedSince{2014_03} */ DEPRECATED(void RemoveSegmentationFromContourList(mitk::Image *segmentation)); /** * @brief Remove interpolation session * @param segmentationImage the session to be removed */ void RemoveInterpolationSession(mitk::Image::Pointer segmentationImage); /** * Replaces the current interpolation session with a new one. All contours form the old * session will be applied to the new session. This only works if the two images have the * geometry * @param oldSession the session which should be replaced * @param newSession the new session which replaces the old one * @return true it the the replacement was successful, false if not (e.g. the image's geometry differs) */ bool ReplaceInterpolationSession(mitk::Image::Pointer oldSession, mitk::Image::Pointer newSession); /** * @brief Removes all sessions */ void RemoveAllInterpolationSessions(); mitk::Image *GetImage(); /** * @brief Get the Contours at a certain timeStep and layerID. * * @param timeStep Time Step from which to get the contours. * @param layerID Layer from which to get the contours. * @return std::vector Returns contours. */ - ContourPositionInformationList& GetContours(unsigned int timeStep, unsigned int layerID); + ContourPositionInformationList* GetContours(unsigned int timeStep, unsigned int layerID); /** * @brief Trigerred with the "Reinit Interpolation" action. The contours are used to repopulate the * surfaceInterpolator data structures so that interpolation can be performed after reloading data. * * @param contourList List of contours extracted * @param contourPlanes List of planes at which the contours were extracted */ void CompleteReinitialization(const std::vector& contourList, std::vector& contourPlanes); /** * @brief Removes contours of a particular label, at a given time step and layerID. * * @param label Label of contour to remove. * @param timeStep Time step in which to remove the contours. * @param layerID Layer in which the contour should be removed. */ void RemoveContours(mitk::Label::PixelType label, unsigned int timeStep, unsigned int layerID); /** * Estimates the memory which is needed to build up the equationsystem for the interpolation. * \returns The percentage of the real memory which will be used by the interpolation */ double EstimatePortionOfNeededMemory(); /** * Adds Contours from the active Label to the interpolation pipeline */ void AddActiveLabelContoursForInterpolation(mitk::Label::PixelType activeLabel); unsigned int GetNumberOfInterpolationSessions(); /** * @brief Removes the contour for a given plane for the current selected segmenation * @param contourInfo the contour which should be removed * @return true if a contour was found and removed, false if no contour was found */ bool RemoveContour(ContourPositionInformation contourInfo); /** * @brief Get the Segmentation Image Node object * * @return DataNode* returns the DataNode containing the segmentation image. */ mitk::DataNode* GetSegmentationImageNode(); protected: SurfaceInterpolationController(); ~SurfaceInterpolationController() override; template void GetImageBase(itk::Image *input, itk::ImageBase<3>::Pointer &result); private: /** * @brief * * @param caller * @param event */ void OnSegmentationDeleted(const itk::Object *caller, const itk::EventObject &event); /** * @brief Function that removes contours of a particular label when the "Remove Label" event is trigerred in the labelSetImage. * */ void OnRemoveLabel(mitk::Label::PixelType removedLabelValue); /** * @brief When a new contour is added to the pipeline or an existing contour is replaced, * the plane geometry information of that contour is added as a child node to the * current node of the segmentation image. This is useful in the retrieval of contour information * when data is reloaded after saving. * * @param contourInfo contourInfo struct to add to data storage. */ void AddPlaneGeometryNodeToDataStorage(const ContourPositionInformation& contourInfo); /** * @brief Function that toggles active label, when the active label is changed. * */ void OnActiveLabel(mitk::Label::PixelType); /** * @brief Clears the interpolation data structures. Called from CompleteReinitialization(). * */ void ClearInterpolationSession(); /** * @brief Add contour to the interpolation pipeline * * @param contourInfo Contour information to be added * @param reinitializationAction If the contour is coming from a reinitialization process or not */ void AddToInterpolationPipeline(ContourPositionInformation& contourInfo, bool reinitializationAction = false); /** * @brief Function to respond to layer changed * */ void OnLayerChanged(); itk::SmartPointer m_ReduceFilter; itk::SmartPointer m_NormalsFilter; itk::SmartPointer m_InterpolateSurfaceFilter; mitk::Surface::Pointer m_Contours; double m_DistanceImageSpacing; vtkSmartPointer m_PolyData; mitk::DataStorage::Pointer m_DataStorage; ContourContainer m_ListOfInterpolationSessions; ContourListMap m_ListOfContours; mitk::Surface::Pointer m_InterpolationResult; unsigned int m_CurrentNumberOfReducedContours; unsigned int m_NumberOfConnectionsAdded; mitk::Image *m_SelectedSegmentation; std::map m_SegmentationObserverTags; mitk::TimePointType m_CurrentTimePoint; unsigned int m_ContourIndex; unsigned int m_ContourPosIndex; unsigned int m_NumberOfLayersInCurrentSegmentation; mitk::Label::PixelType m_PreviousActiveLabelValue; mitk::Label::PixelType m_CurrentActiveLabelValue; unsigned int m_PreviousLayerIndex; unsigned int m_CurrentLayerIndex; }; namespace ContourExt { /** * @brief Returns the plane the contour belongs to. * * @param ContourNormal * @return size_t */ size_t GetContourOrientation(const mitk::Vector3D& ContourNormal); /** * @brief Function used to compute an interior point of the contour. * Used to react to the merge label and erase label actions. * * * @tparam VImageDimension Dimension of the image * @param contour Contour for which to compute the interior point * @param labelSetImage Label Set Image For which to find the contour * @param currentTimePoint Current Time Point of the Image * @return mitk::Point3D The returned point in the interior of the contour.s */ template mitk::Point3D ComputeInteriorPointOfContour(const mitk::SurfaceInterpolationController::ContourPositionInformation& contour, mitk::LabelSetImage * labelSetImage, mitk::TimePointType currentTimePoint); /** * @brief Get a Grid points within the bounding box of the contour at a certain spacing. * * @param planeDimension Plane orientation (Sagittal, Coronal, Axial) * @param startDim1 Starting coordinate along dimension 1 to start the grid point sampling from * @param numPointsToSampleDim1 Number of points to sample along dimension 1 * @param deltaDim1 Spacing for dimension 1 at which points should be sampled * @param startDim2 Starting coordinate along dimension 2 to start the grid point sampling from * @param numPointsToSampleDim2 Number of points to sample along dimension 2 * @param deltaDim2 Spacing for dimension 1 at which points should be sampled * @param valuePlaneDim Slice index of the plane in the volume * @return std::vector< mitk::Point3D > The computed grid points are returned by the function. */ std::vector< mitk::Point3D > GetBoundingBoxGridPoints(size_t planeDimension, double startDim1, size_t numPointsToSampleDim1, double deltaDim1, double startDim2, size_t numPointsToSampleDim2, double deltaDim2, double valuePlaneDim); }; } #endif