diff --git a/Modules/CEST/files.cmake b/Modules/CEST/files.cmake index 4655b5010f..19a794674c 100644 --- a/Modules/CEST/files.cmake +++ b/Modules/CEST/files.cmake @@ -1,10 +1,11 @@ file(GLOB_RECURSE H_FILES RELATIVE "${CMAKE_CURRENT_SOURCE_DIR}" "${CMAKE_CURRENT_SOURCE_DIR}/include/*") set(CPP_FILES mitkCESTImageNormalizationFilter.cpp mitkCustomTagParser.cpp ) set(RESOURCE_FILES 1416.json + 1485.json ) diff --git a/Modules/CEST/include/mitkCustomTagParser.h b/Modules/CEST/include/mitkCustomTagParser.h index ef5671e2a0..7a53928ef1 100644 --- a/Modules/CEST/include/mitkCustomTagParser.h +++ b/Modules/CEST/include/mitkCustomTagParser.h @@ -1,120 +1,123 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #ifndef MITKCUSTOMTAGPARSER_H #define MITKCUSTOMTAGPARSER_H #include #include #include namespace mitk { /** The custom tag parser can be used to parse the custom dicom tag of the siemens private tag (0x0029, 0x1020) to extract relevant CEST data. An initial parsing determines whether the provided string belongs to CEST data at all. If the "tSequenceFileName" is of the format "{WHATEVER}CEST_Rev####" it is assumed that the data is indeed CEST data and was taken with revision #### (not limited to four digits). Which custom parameters to save and to which property name can be controlled by a json file. This file can be either provided as a resource for the MitkCEST module during compilation or placed next to the MitkCEST library in your binary folder. The expected format for the file "REVISIONNUMBER.json":
{
"REVISIONNUMBER" : "revision_json",
"sWiPMemBlock.alFree[1]" : "AdvancedMode",
"sWiPMemBlock.alFree[2]" : "RetreatMode"
}
where :
  • REVISIONNUMBER is the revision number of this json parameter mapping (files with non digit characters in their name will be ignored)
  • sWiPMemBlock.alFree[1] is the name of one parameter in the private dicom tag
  • AdvancedMode is the name of the property the content of sWiPMemBlock.alFree[1] should be saved to
\note It is assumed that the entire content of tag (0x0029, 0x1020) is provided and that it es hex encoded (12\23\04...). If the sampling type is list it will try to access LIST.txt at the location provided in the constructor to read the offsets. */ class MITKCEST_EXPORT CustomTagParser { public: /// the constructor expects a path to one of the files to be loaded or the directory of the dicom files CustomTagParser(std::string relevantFile); /// parse the provided dicom property and return a property list based on the closest revision parameter mapping mitk::PropertyList::Pointer ParseDicomProperty(mitk::TemporoSpatialStringProperty *dicomProperty); /// parse the provided string and return a property list based on the closest revision parameter mapping mitk::PropertyList::Pointer ParseDicomPropertyString(std::string dicomPropertyString); void SetParseStrategy(std::string parseStrategy); void SetRevisionMappingStrategy(std::string revisionMappingStrategy); /// name of the property for the offsets, including normalization offsets static const std::string m_OffsetsPropertyName; /// name of the property for the data acquisition revision static const std::string m_RevisionPropertyName; /// name of the property for the json parameter mapping revision static const std::string m_JSONRevisionPropertyName; /// prefix for all CEST related property names static const std::string m_CESTPropertyPrefix; protected: std::string GetRevisionAppropriateJSONString(std::string revisionString); void GetClosestLowerRevision(std::string revisionString); std::string GetClosestLowerRevision(std::string revisionString, std::vector availableRevisionsVector); + + /// Decides whether or not the image is likely to be a T1Map, if not it is assumed to be a CEST sequence + bool IsT1Sequence(std::string preparationType, std::string recoveryMode, std::string spoilingType, std::string revisionString); /// Get a string filled with the properly formated offsets based on the sampling type and offset std::string GetOffsetString(std::string samplingType, std::string offset, std::string measurements); /// returns a vector revision numbers of all REVISIONNUMBER.json found beside the MitkCEST library std::vector GetExternalRevisions(); /// returns a vector revision numbers of all REVISIONNUMBER.json provided as resources during the compile std::vector GetInternalRevisions(); /// returns the path where external jsons are expected to be located std::string GetExternalJSONDirectory(); /// the closest lower revision provided as resource, empty if none found std::string m_ClosestInternalRevision; /// the closest lower revision provided as a json beside the library, empty if none found std::string m_ClosestExternalRevision; /// revision independent mapping to inject into the revision dependent json string static const std::string m_RevisionIndependentMapping; /// default revision dependent json string if none is found static const std::string m_DefaultJsonString; /// path to the dicom data std::string m_DicomDataPath; /// Should the kind of data be automatically determined or should it be parsed as a specific one std::string m_ParseStrategy; /// How to handle parameter mapping based on absent revision jsons std::string m_RevisionMappingStrategy; }; } #endif // MITKCUSTOMTAGPARSER_H diff --git a/Modules/CEST/resource/1485.json b/Modules/CEST/resource/1485.json new file mode 100644 index 0000000000..3119952174 --- /dev/null +++ b/Modules/CEST/resource/1485.json @@ -0,0 +1,24 @@ +{ + "1485" : "revision_json", + "sWiPMemBlock.alFree[1]" : "AdvancedMode", + "sWiPMemBlock.alFree[2]" : "RecoveryMode", + "sWiPMemBlock.alFree[3]" : "DoubleIrrMode", + "sWiPMemBlock.alFree[4]" : "MtMode", + "sWiPMemBlock.alFree[5]" : "PreparationType", + "sWiPMemBlock.alFree[6]" : "PulseType", + "sWiPMemBlock.alFree[7]" : "SamplingType", + "sWiPMemBlock.alFree[8]" : "SpoilingType", + "sWiPMemBlock.alFree[9]" : "measurements", + "sWiPMemBlock.alFree[10]" : "NumberRFBlocks", + "sWiPMemBlock.alFree[11]" : "PulsesPerRFBlock", + "sWiPMemBlock.alFree[12]" : "PulseDuration", + "sWiPMemBlock.alFree[13]" : "DutyCycle", + "sWiPMemBlock.alFree[14]" : "RecoveryTime", + "sWiPMemBlock.alFree[15]" : "RecoveryTimeM0", + "sWiPMemBlock.adFree[1]" : "Offset", + "sWiPMemBlock.adFree[2]" : "B1Amplitude", + "sWiPMemBlock.adFree[3]" : "AdiabaticPulseMu", + "sWiPMemBlock.adFree[4]" : "AdiabaticPulseBW", + "sWiPMemBlock.adFree[5]" : "AdiabaticPulseLength", + "sWiPMemBlock.adFree[6]" : "AdiabaticPulseAmp" +} diff --git a/Modules/CEST/test/files.cmake b/Modules/CEST/test/files.cmake index 88653735a2..46c8a5ee61 100644 --- a/Modules/CEST/test/files.cmake +++ b/Modules/CEST/test/files.cmake @@ -1,6 +1,7 @@ set(MODULE_TESTS mitkCustomTagParserTest.cpp + mitkCESTDICOMReaderServiceTest.cpp ) SET(MODULE_CUSTOM_TESTS ) diff --git a/Modules/CEST/test/mitkCESTDICOMReaderServiceTest.cpp b/Modules/CEST/test/mitkCESTDICOMReaderServiceTest.cpp new file mode 100644 index 0000000000..3f0f946877 --- /dev/null +++ b/Modules/CEST/test/mitkCESTDICOMReaderServiceTest.cpp @@ -0,0 +1,81 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +// Testing +#include "mitkTestFixture.h" +#include "mitkTestingMacros.h" + +// std includes +#include + +// MITK includes +#include +#include +#include + + +// VTK includes +#include + +class mitkCESTDICOMReaderServiceTestSuite : public mitk::TestFixture +{ + CPPUNIT_TEST_SUITE(mitkCESTDICOMReaderServiceTestSuite); + + // Test the dicom property parsing + MITK_TEST(LoadCESTDICOMData_Success); + MITK_TEST(LoadT1DICOMData_Success); + + CPPUNIT_TEST_SUITE_END(); + +private: + + +public: + void setUp() override + { + + } + + void tearDown() override + { + } + + void LoadCESTDICOMData_Success() + { + mitk::IFileReader::Options options; + options["Force type"] = std::string( "Automatic" ); + options["Revision mapping"] = std::string( "Strict" ); + + mitk::Image::Pointer cestImage = mitk::IOUtil::Load(GetTestDataFilePath("CEST/B1=0.6MUT/MEI_NER_PHANTOM.MR.E0202_MEISSNER.0587.0001.2017.10.25.22.11.10.373351.41828677.IMA"), options); + CPPUNIT_ASSERT_MESSAGE("Make certain offsets have been correctly loaded for CEST image." ,cestImage->GetProperty("CEST.Offsets")->GetValueAsString() == "-300 2 -2 1.92982 -1.92982 1.85965 -1.85965 1.78947 -1.78947 1.7193 -1.7193 1.64912 -1.64912 1.57895 -1.57895 1.50877 -1.50877 1.4386 -1.4386 1.36842 -1.36842 1.29825 -1.29825 1.22807 -1.22807 1.15789 -1.15789 1.08772 -1.08772 1.01754 -1.01754 0.947368 -0.947368 0.877193 -0.877193 0.807018 -0.807018 0.736842 -0.736842 0.666667 -0.666667 0.596491 -0.596491 0.526316 -0.526316 0.45614 -0.45614 0.385965 -0.385965 0.315789 -0.315789 0.245614 -0.245614 0.175439 -0.175439 0.105263 -0.105263 0.0350877 -0.0350877"); + std::string temp; + CPPUNIT_ASSERT_MESSAGE("Make certain image is not loaded as T1.", !cestImage->GetPropertyList()->GetStringProperty("CEST.TREC", temp)); + } + + void LoadT1DICOMData_Success() + { + mitk::IFileReader::Options options; + options["Force type"] = std::string("Automatic"); + options["Revision mapping"] = std::string("Strict"); + + mitk::Image::Pointer cestImage = mitk::IOUtil::Load(GetTestDataFilePath("CEST/T1MAP/MEI_NER_PHANTOM.MR.E0202_MEISSNER.0279.0001.2017.10.25.20.21.27.996834.41803047.IMA"), options); + std::string temp; + CPPUNIT_ASSERT_MESSAGE("Make certain image is loaded as T1.", cestImage->GetPropertyList()->GetStringProperty("CEST.TREC", temp)); + } + +}; + +MITK_TEST_SUITE_REGISTRATION(mitkCESTDICOMReaderService) diff --git a/Plugins/org.mitk.gui.qt.cest/src/internal/QmitkCESTStatisticsView.cpp b/Plugins/org.mitk.gui.qt.cest/src/internal/QmitkCESTStatisticsView.cpp index 3c4f8f51fb..4040af43a7 100644 --- a/Plugins/org.mitk.gui.qt.cest/src/internal/QmitkCESTStatisticsView.cpp +++ b/Plugins/org.mitk.gui.qt.cest/src/internal/QmitkCESTStatisticsView.cpp @@ -1,899 +1,956 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ //itk #include "itksys/SystemTools.hxx" #include #include // Blueberry #include #include // Qmitk #include "QmitkCESTStatisticsView.h" // Qt #include #include // qwt #include // mitk #include #include #include #include #include #include #include #include #include #include // boost #include #include //stl #include #include #include #include #include #include +namespace +{ + template + void GetSortPermutation(std::vector &out, + const std::vector &determiningVector, + Compare compare = std::less()) + { + out.resize(determiningVector.size()); + std::iota(out.begin(), out.end(), 0); + + std::sort(out.begin(), out.end(), [&](unsigned i, unsigned j) { + return compare(determiningVector[i], determiningVector[j]); + }); + } + + template + void ApplyPermutation(const std::vector &order, std::vector &vectorToSort) + { + assert(order.size() == vectorToSort.size()); + std::vector tempVector(vectorToSort.size()); + for (unsigned i = 0; i < vectorToSort.size(); i++) + { + tempVector[i] = vectorToSort[order[i]]; + } + vectorToSort = tempVector; + } + + template + void ApplyPermutation(const std::vector &order, + std::vector ¤tVector, + std::vector &... remainingVectors) + { + ApplyPermutation(order, currentVector); + ApplyPermutation(order, remainingVectors...); + } + + template + void SortVectors(const std::vector &orderDeterminingVector, + Compare comparison, + std::vector &... vectorsToBeSorted) + { + std::vector order; + GetSortPermutation(order, orderDeterminingVector, comparison); + ApplyPermutation(order, vectorsToBeSorted...); + } +} + const std::string QmitkCESTStatisticsView::VIEW_ID = "org.mitk.views.ceststatistics"; static const int STAT_TABLE_BASE_HEIGHT = 180; QmitkCESTStatisticsView::QmitkCESTStatisticsView(QObject* /*parent*/, const char* /*name*/) { this->m_CalculatorThread = new QmitkImageStatisticsCalculationThread; m_currentSelectedPosition.Fill(0.0); m_currentSelectedTimeStep = 0; m_CrosshairPointSet = mitk::PointSet::New(); } QmitkCESTStatisticsView::~QmitkCESTStatisticsView() { while (this->m_CalculatorThread->isRunning()) // wait until thread has finished { itksys::SystemTools::Delay(100); } delete this->m_CalculatorThread; } void QmitkCESTStatisticsView::SetFocus() { m_Controls.threeDimToFourDimPushButton->setFocus(); } void QmitkCESTStatisticsView::CreateQtPartControl( QWidget *parent ) { // create GUI widgets from the Qt Designer's .ui file m_Controls.setupUi( parent ); connect(m_Controls.threeDimToFourDimPushButton, SIGNAL(clicked()), this, SLOT(OnThreeDimToFourDimPushButtonClicked())); connect((QObject*) this->m_CalculatorThread, SIGNAL(finished()), this, SLOT(OnThreadedStatisticsCalculationEnds()), Qt::QueuedConnection); connect((QObject*)(this->m_Controls.m_CopyStatisticsToClipboardPushButton), SIGNAL(clicked()), (QObject*) this, SLOT(OnCopyStatisticsToClipboardPushButtonClicked())); connect((QObject*)(this->m_Controls.normalizeImagePushButton), SIGNAL(clicked()), (QObject*) this, SLOT(OnNormalizeImagePushButtonClicked())); connect((QObject*)(this->m_Controls.fixedRangeCheckBox), SIGNAL(toggled(bool)), (QObject*) this, SLOT(OnFixedRangeCheckBoxToggled(bool))); connect((QObject*)(this->m_Controls.fixedRangeLowerDoubleSpinBox), SIGNAL(editingFinished()), (QObject*) this, SLOT(OnFixedRangeDoubleSpinBoxChanged())); connect((QObject*)(this->m_Controls.fixedRangeUpperDoubleSpinBox), SIGNAL(editingFinished()), (QObject*) this, SLOT(OnFixedRangeDoubleSpinBoxChanged())); m_Controls.normalizeImagePushButton->setEnabled(false); m_Controls.threeDimToFourDimPushButton->setEnabled(false); this->m_SliceChangeListener.RenderWindowPartActivated(this->GetRenderWindowPart()); connect(&m_SliceChangeListener, SIGNAL(SliceChanged()), this, SLOT(OnSliceChanged())); } void QmitkCESTStatisticsView::RenderWindowPartActivated(mitk::IRenderWindowPart* renderWindowPart) { this->m_SliceChangeListener.RenderWindowPartActivated(renderWindowPart); } void QmitkCESTStatisticsView::RenderWindowPartDeactivated( mitk::IRenderWindowPart* renderWindowPart) { this->m_SliceChangeListener.RenderWindowPartDeactivated(renderWindowPart); } void QmitkCESTStatisticsView::OnSelectionChanged( berry::IWorkbenchPart::Pointer /*source*/, const QList& nodes ) { if (nodes.empty()) { std::stringstream message; message << "Please select an image."; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); this->Clear(); return; } // iterate all selected objects bool atLeastOneWasCESTImage = false; foreach( mitk::DataNode::Pointer node, nodes ) { if (node.IsNull()) { continue; } if( dynamic_cast(node->GetData()) != nullptr ) { m_Controls.labelWarning->setVisible( false ); bool zSpectrumSet = SetZSpectrum(dynamic_cast(node->GetData()->GetProperty(mitk::CustomTagParser::m_OffsetsPropertyName.c_str()).GetPointer())); atLeastOneWasCESTImage = atLeastOneWasCESTImage || zSpectrumSet; if (zSpectrumSet) { m_ZImage = dynamic_cast(node->GetData()); } else { m_MaskImage = dynamic_cast(node->GetData()); } } if (dynamic_cast(node->GetData()) != nullptr) { m_MaskPlanarFigure = dynamic_cast(node->GetData()); } if (dynamic_cast(node->GetData()) != nullptr) { m_PointSet = dynamic_cast(node->GetData()); } } // We only want to offer normalization or timestep copying if one object is selected if (nodes.size() == 1) { if (dynamic_cast(nodes.front()->GetData()) ) { m_Controls.normalizeImagePushButton->setEnabled(atLeastOneWasCESTImage); m_Controls.threeDimToFourDimPushButton->setDisabled(atLeastOneWasCESTImage); } else { m_Controls.normalizeImagePushButton->setEnabled(false); m_Controls.threeDimToFourDimPushButton->setEnabled(false); std::stringstream message; message << "The selected node is not an image."; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); } this->Clear(); return; } // we always need a mask, either image or planar figure as well as an image for further processing if (nodes.size() != 2) { this->Clear(); return; } m_Controls.normalizeImagePushButton->setEnabled(false); m_Controls.threeDimToFourDimPushButton->setEnabled(false); if (!atLeastOneWasCESTImage) { std::stringstream message; message << "None of the selected data nodes contains required CEST meta information"; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); this->Clear(); return; } bool bothAreImages = (m_ZImage.GetPointer() != nullptr) && (m_MaskImage.GetPointer() != nullptr); if (bothAreImages) { bool geometriesMatch = mitk::Equal(*(m_ZImage->GetTimeGeometry()), *(m_MaskImage->GetTimeGeometry()), mitk::eps, false); if (!geometriesMatch) { std::stringstream message; message << "The selected images have different geometries."; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); this->Clear(); return; } } if (!this->DataSanityCheck()) { this->Clear(); return; } if (m_PointSet.IsNull()) { // initialize thread and trigger it this->m_CalculatorThread->SetIgnoreZeroValueVoxel(false); this->m_CalculatorThread->Initialize(m_ZImage, m_MaskImage, m_MaskPlanarFigure); std::stringstream message; message << "Calculating statistics..."; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); try { // Compute statistics this->m_CalculatorThread->start(); } catch (const mitk::Exception& e) { std::stringstream message; message << "" << e.GetDescription() << ""; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); } catch (const std::runtime_error &e) { // In case of exception, print error message on GUI std::stringstream message; message << "" << e.what() << ""; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); } catch (const std::exception &e) { MITK_ERROR << "Caught exception: " << e.what(); // In case of exception, print error message on GUI std::stringstream message; message << "Error! Unequal Dimensions of Image and Segmentation. No recompute possible "; m_Controls.labelWarning->setText(message.str().c_str()); m_Controls.labelWarning->show(); } while (this->m_CalculatorThread->isRunning()) // wait until thread has finished { itksys::SystemTools::Delay(100); } } if (m_PointSet.IsNotNull()) { if (m_ZImage->GetDimension() == 4) { AccessFixedDimensionByItk(m_ZImage, PlotPointSet, 4); } else { MITK_WARN << "Expecting a 4D image."; } } } void QmitkCESTStatisticsView::OnThreadedStatisticsCalculationEnds() { this->m_Controls.m_DataViewWidget->SetAxisTitle(QwtPlot::Axis::xBottom, "delta w"); this->m_Controls.m_DataViewWidget->SetAxisTitle(QwtPlot::Axis::yLeft, "z"); const std::vector &statistics = this->m_CalculatorThread->GetStatisticsData(); QmitkPlotWidget::DataVector::size_type numberOfSpectra = this->m_zSpectrum.size(); QmitkPlotWidget::DataVector means(numberOfSpectra); QmitkPlotWidget::DataVector stdevs(numberOfSpectra); for (unsigned int index = 0; index < numberOfSpectra; ++index) { means[index] = statistics[index]->GetMean(); stdevs[index] = statistics[index]->GetStd(); } QmitkPlotWidget::DataVector xValues = this->m_zSpectrum; RemoveMZeros(xValues, means, stdevs); + ::SortVectors(xValues, std::less(), xValues, means, stdevs); unsigned int curveId = this->m_Controls.m_DataViewWidget->InsertCurve("Spectrum"); this->m_Controls.m_DataViewWidget->SetCurveData(curveId, xValues, means, stdevs, stdevs); this->m_Controls.m_DataViewWidget->SetErrorPen(curveId, QPen(Qt::blue)); QwtSymbol* blueSymbol = new QwtSymbol(QwtSymbol::Rect, QColor(Qt::blue), QColor(Qt::blue), QSize(8, 8)); this->m_Controls.m_DataViewWidget->SetCurveSymbol(curveId, blueSymbol); this->m_Controls.m_DataViewWidget->SetLegendAttribute(curveId, QwtPlotCurve::LegendShowSymbol); QwtLegend* legend = new QwtLegend(); legend->setFrameShape(QFrame::Box); legend->setFrameShadow(QFrame::Sunken); legend->setLineWidth(1); this->m_Controls.m_DataViewWidget->SetLegend(legend, QwtPlot::BottomLegend); m_Controls.m_DataViewWidget->GetPlot()->axisScaleEngine(QwtPlot::Axis::xBottom)->setAttributes(QwtScaleEngine::Inverted); this->m_Controls.m_DataViewWidget->Replot(); m_Controls.labelWarning->setVisible(false); if (this->m_Controls.fixedRangeCheckBox->isChecked()) { this->m_Controls.m_DataViewWidget->GetPlot()->setAxisAutoScale(2, false); this->m_Controls.m_DataViewWidget->GetPlot()->setAxisScale(2, this->m_Controls.fixedRangeLowerDoubleSpinBox->value(), this->m_Controls.fixedRangeUpperDoubleSpinBox->value()); } else { this->m_Controls.m_DataViewWidget->GetPlot()->setAxisAutoScale(2, true); } if(this->DataSanityCheck()) { this->FillStatisticsTableView(this->m_CalculatorThread->GetStatisticsData(), this->m_CalculatorThread->GetStatisticsImage()); } else { this->Clear(); } } void QmitkCESTStatisticsView::OnFixedRangeDoubleSpinBoxChanged() { if (this->m_Controls.fixedRangeCheckBox->isChecked()) { this->m_Controls.m_DataViewWidget->GetPlot()->setAxisAutoScale(2, false); this->m_Controls.m_DataViewWidget->GetPlot()->setAxisScale(2, this->m_Controls.fixedRangeLowerDoubleSpinBox->value(), this->m_Controls.fixedRangeUpperDoubleSpinBox->value()); } this->m_Controls.m_DataViewWidget->Replot(); } template void QmitkCESTStatisticsView::PlotPointSet(itk::Image* image) { this->m_Controls.m_DataViewWidget->SetAxisTitle(QwtPlot::Axis::xBottom, "delta w"); this->m_Controls.m_DataViewWidget->SetAxisTitle(QwtPlot::Axis::yLeft, "z"); QmitkPlotWidget::DataVector::size_type numberOfSpectra = this->m_zSpectrum.size(); mitk::PointSet::Pointer internalPointset; if (m_PointSet.IsNotNull()) { internalPointset = m_PointSet; } else { internalPointset = m_CrosshairPointSet; } if (internalPointset.IsNull()) { return; } + if (!this->DataSanityCheck()) + { + m_Controls.labelWarning->setText("Data can not be plotted, internally inconsistent."); + m_Controls.labelWarning->show(); + return; + } + auto maxIndex = internalPointset->GetMaxId().Index(); for (std::size_t number = 0; number < maxIndex + 1; ++number) { mitk::PointSet::PointType point; if (!internalPointset->GetPointIfExists(number, &point)) { continue; } if (!this->m_ZImage->GetGeometry()->IsInside(point)) { continue; } itk::Index<3> itkIndex; this->m_ZImage->GetGeometry()->WorldToIndex(point, itkIndex); itk::Index itkIndexTime; itkIndexTime[0] = itkIndex[0]; itkIndexTime[1] = itkIndex[1]; itkIndexTime[2] = itkIndex[2]; QmitkPlotWidget::DataVector values(numberOfSpectra); for (std::size_t step = 0; step < numberOfSpectra; ++step) { if( VImageDimension == 4 ) { itkIndexTime[3] = step; } values[step] = image->GetPixel(itkIndexTime); } std::stringstream name; name << "Point " << number; // Qcolor enums go from 0 to 19, but 19 is transparent and 0,1 are for bitmaps // 3 is white and thus not visible QColor color(static_cast(number % 17 + 4)); QmitkPlotWidget::DataVector xValues = this->m_zSpectrum; RemoveMZeros(xValues, values); + ::SortVectors(xValues, std::less(), xValues, values); unsigned int curveId = this->m_Controls.m_DataViewWidget->InsertCurve(name.str().c_str()); this->m_Controls.m_DataViewWidget->SetCurveData(curveId, xValues, values); this->m_Controls.m_DataViewWidget->SetCurvePen(curveId, QPen(color)); QwtSymbol* symbol = new QwtSymbol(QwtSymbol::Rect, color, color, QSize(8, 8)); this->m_Controls.m_DataViewWidget->SetCurveSymbol(curveId, symbol); this->m_Controls.m_DataViewWidget->SetLegendAttribute(curveId, QwtPlotCurve::LegendShowSymbol); } if (this->m_Controls.fixedRangeCheckBox->isChecked()) { this->m_Controls.m_DataViewWidget->GetPlot()->setAxisAutoScale(2, false); this->m_Controls.m_DataViewWidget->GetPlot()->setAxisScale(2, this->m_Controls.fixedRangeLowerDoubleSpinBox->value(), this->m_Controls.fixedRangeUpperDoubleSpinBox->value()); } else { this->m_Controls.m_DataViewWidget->GetPlot()->setAxisAutoScale(2, true); } QwtLegend* legend = new QwtLegend(); legend->setFrameShape(QFrame::Box); legend->setFrameShadow(QFrame::Sunken); legend->setLineWidth(1); this->m_Controls.m_DataViewWidget->SetLegend(legend, QwtPlot::BottomLegend); m_Controls.m_DataViewWidget->GetPlot()->axisScaleEngine(QwtPlot::Axis::xBottom)->setAttributes(QwtScaleEngine::Inverted); this->m_Controls.m_DataViewWidget->Replot(); m_Controls.labelWarning->setVisible(false); } void QmitkCESTStatisticsView::OnFixedRangeCheckBoxToggled(bool state) { this->m_Controls.fixedRangeLowerDoubleSpinBox->setEnabled(state); this->m_Controls.fixedRangeUpperDoubleSpinBox->setEnabled(state); } void QmitkCESTStatisticsView::OnNormalizeImagePushButtonClicked() { QList nodes = this->GetDataManagerSelection(); if (nodes.empty()) return; mitk::DataNode* node = nodes.front(); if (!node) { // Nothing selected. Inform the user and return QMessageBox::information(nullptr, "CEST View", "Please load and select an image before starting image processing."); return; } // here we have a valid mitk::DataNode // a node itself is not very useful, we need its data item (the image) mitk::BaseData* data = node->GetData(); if (data) { // test if this data item is an image or not (could also be a surface or something totally different) mitk::Image* image = dynamic_cast(data); if (image) { std::string offsets = ""; bool hasOffsets = image->GetPropertyList()->GetStringProperty( mitk::CustomTagParser::m_OffsetsPropertyName.c_str() ,offsets); if (!hasOffsets) { QMessageBox::information(nullptr, "CEST View", "Selected image was missing CEST offset information."); return; } if (image->GetDimension() == 4) { auto normalizationFilter = mitk::CESTImageNormalizationFilter::New(); normalizationFilter->SetInput(image); normalizationFilter->Update(); auto resultImage = normalizationFilter->GetOutput(); mitk::DataNode::Pointer dataNode = mitk::DataNode::New(); dataNode->SetData(resultImage); std::string normalizedName = node->GetName() + "_normalized"; dataNode->SetName(normalizedName); this->GetDataStorage()->Add(dataNode); } this->Clear(); } } } void QmitkCESTStatisticsView::RemoveMZeros(QmitkPlotWidget::DataVector& xValues, QmitkPlotWidget::DataVector& yValues) { QmitkPlotWidget::DataVector tempX; QmitkPlotWidget::DataVector tempY; for (std::size_t index = 0; index < xValues.size(); ++index) { if ((xValues.at(index) < -299) || (xValues.at(index)) > 299) { // do not include } else { tempX.push_back(xValues.at(index)); tempY.push_back(yValues.at(index)); } } xValues = tempX; yValues = tempY; } void QmitkCESTStatisticsView::RemoveMZeros(QmitkPlotWidget::DataVector& xValues, QmitkPlotWidget::DataVector& yValues, QmitkPlotWidget::DataVector& stdDevs) { QmitkPlotWidget::DataVector tempX; QmitkPlotWidget::DataVector tempY; QmitkPlotWidget::DataVector tempDevs; for (std::size_t index = 0; index < xValues.size(); ++index) { if ((xValues.at(index) < -299) || (xValues.at(index)) > 299) { // do not include } else { tempX.push_back(xValues.at(index)); tempY.push_back(yValues.at(index)); tempDevs.push_back(stdDevs.at(index)); } } xValues = tempX; yValues = tempY; stdDevs = tempDevs; } void QmitkCESTStatisticsView::OnThreeDimToFourDimPushButtonClicked() { QList nodes = this->GetDataManagerSelection(); if (nodes.empty()) return; mitk::DataNode* node = nodes.front(); if (!node) { // Nothing selected. Inform the user and return QMessageBox::information( nullptr, "CEST View", "Please load and select an image before starting image processing."); return; } // here we have a valid mitk::DataNode // a node itself is not very useful, we need its data item (the image) mitk::BaseData* data = node->GetData(); if (data) { // test if this data item is an image or not (could also be a surface or something totally different) mitk::Image* image = dynamic_cast( data ); if (image) { if (image->GetDimension() == 4) { AccessFixedDimensionByItk(image, CopyTimesteps, 4); } this->Clear(); } } } template void QmitkCESTStatisticsView::CopyTimesteps(itk::Image* image) { typedef itk::Image ImageType; //typedef itk::PasteImageFilter PasteImageFilterType; unsigned int numberOfTimesteps = image->GetLargestPossibleRegion().GetSize(3); typename ImageType::RegionType sourceRegion = image->GetLargestPossibleRegion(); sourceRegion.SetSize(3, 1); typename ImageType::RegionType targetRegion = image->GetLargestPossibleRegion(); targetRegion.SetSize(3, 1); for (unsigned int timestep = 1; timestep < numberOfTimesteps; ++timestep) { targetRegion.SetIndex(3, timestep); itk::ImageRegionConstIterator sourceIterator(image, sourceRegion); itk::ImageRegionIterator targetIterator(image, targetRegion); while (!sourceIterator.IsAtEnd()) { targetIterator.Set(sourceIterator.Get()); ++sourceIterator; ++targetIterator; } } } bool QmitkCESTStatisticsView::SetZSpectrum(mitk::StringProperty* zSpectrumProperty) { if (zSpectrumProperty == nullptr) { return false; } mitk::LocaleSwitch localeSwitch("C"); std::string zSpectrumString = zSpectrumProperty->GetValueAsString(); std::istringstream iss(zSpectrumString); std::vector zSpectra; std::copy(std::istream_iterator(iss), std::istream_iterator(), std::back_inserter(zSpectra)); m_zSpectrum.clear(); m_zSpectrum.resize(0); for (auto const &spectrumString : zSpectra) { m_zSpectrum.push_back(std::stod(spectrumString)); } return (m_zSpectrum.size() > 0); } void QmitkCESTStatisticsView::FillStatisticsTableView( const std::vector &s, const mitk::Image *image) { this->m_Controls.m_StatisticsTable->setColumnCount(image->GetTimeSteps()); this->m_Controls.m_StatisticsTable->horizontalHeader()->setVisible(image->GetTimeSteps() > 1); int decimals = 2; mitk::PixelType doublePix = mitk::MakeScalarPixelType< double >(); mitk::PixelType floatPix = mitk::MakeScalarPixelType< float >(); if (image->GetPixelType() == doublePix || image->GetPixelType() == floatPix) { decimals = 5; } for (unsigned int t = 0; t < image->GetTimeSteps(); t++) { this->m_Controls.m_StatisticsTable->setHorizontalHeaderItem(t, new QTableWidgetItem(QString::number(m_zSpectrum[t]))); this->m_Controls.m_StatisticsTable->setItem(0, t, new QTableWidgetItem( QString("%1").arg(s[t]->GetMean(), 0, 'f', decimals))); this->m_Controls.m_StatisticsTable->setItem(1, t, new QTableWidgetItem( QString("%1").arg(s[t]->GetStd(), 0, 'f', decimals))); this->m_Controls.m_StatisticsTable->setItem(2, t, new QTableWidgetItem( QString("%1").arg(s[t]->GetRMS(), 0, 'f', decimals))); QString max; max.append(QString("%1").arg(s[t]->GetMax(), 0, 'f', decimals)); max += " ("; for (unsigned int i = 0; iGetMaxIndex().size(); i++) { max += QString::number(s[t]->GetMaxIndex()[i]); if (iGetMaxIndex().size() - 1) max += ","; } max += ")"; this->m_Controls.m_StatisticsTable->setItem(3, t, new QTableWidgetItem(max)); QString min; min.append(QString("%1").arg(s[t]->GetMin(), 0, 'f', decimals)); min += " ("; for (unsigned int i = 0; iGetMinIndex().size(); i++) { min += QString::number(s[t]->GetMinIndex()[i]); if (iGetMinIndex().size() - 1) min += ","; } min += ")"; this->m_Controls.m_StatisticsTable->setItem(4, t, new QTableWidgetItem(min)); this->m_Controls.m_StatisticsTable->setItem(5, t, new QTableWidgetItem( QString("%1").arg(s[t]->GetN()))); const mitk::BaseGeometry *geometry = image->GetGeometry(); if (geometry != nullptr) { const mitk::Vector3D &spacing = image->GetGeometry()->GetSpacing(); double volume = spacing[0] * spacing[1] * spacing[2] * (double)s[t]->GetN(); this->m_Controls.m_StatisticsTable->setItem(6, t, new QTableWidgetItem( QString("%1").arg(volume, 0, 'f', decimals))); } else { this->m_Controls.m_StatisticsTable->setItem(6, t, new QTableWidgetItem( "NA")); } } this->m_Controls.m_StatisticsTable->resizeColumnsToContents(); int height = STAT_TABLE_BASE_HEIGHT; if (this->m_Controls.m_StatisticsTable->horizontalHeader()->isVisible()) height += this->m_Controls.m_StatisticsTable->horizontalHeader()->height(); //if (this->m_Controls.m_StatisticsTable->horizontalScrollBar()->isVisible()) // height += this->m_Controls.m_StatisticsTable->horizontalScrollBar()->height(); this->m_Controls.m_StatisticsTable->setMinimumHeight(height); this->m_Controls.m_StatisticsGroupBox->setEnabled(true); this->m_Controls.m_StatisticsTable->setEnabled(true); } void QmitkCESTStatisticsView::InvalidateStatisticsTableView() { this->m_Controls.m_StatisticsTable->horizontalHeader()->setVisible(false); this->m_Controls.m_StatisticsTable->setColumnCount(1); for (int i = 0; i < this->m_Controls.m_StatisticsTable->rowCount(); ++i) { { this->m_Controls.m_StatisticsTable->setItem(i, 0, new QTableWidgetItem("NA")); } } this->m_Controls.m_StatisticsTable->setMinimumHeight(STAT_TABLE_BASE_HEIGHT); this->m_Controls.m_StatisticsTable->setEnabled(false); } bool QmitkCESTStatisticsView::DataSanityCheck() { QmitkPlotWidget::DataVector::size_type numberOfSpectra = m_zSpectrum.size(); // if we do not have a spectrum, the data can not be processed if (numberOfSpectra == 0) { return false; } // if we do not have CEST image data, the data can not be processed if (m_ZImage.IsNull()) { return false; } // if the CEST image data and the meta information do not match, the data can not be processed if (numberOfSpectra != m_ZImage->GetTimeSteps()) { + MITK_INFO << "CEST meta information and number of volumes does not match."; return false; } // if we have neither a mask image, a point set nor a mask planar figure, we can not do statistics // statistics on the whole image would not make sense - if (m_MaskImage.IsNull() && m_MaskPlanarFigure.IsNull() && m_PointSet.IsNull() ) + if (m_MaskImage.IsNull() && m_MaskPlanarFigure.IsNull() && m_PointSet.IsNull() && m_CrosshairPointSet->IsEmpty()) { return false; } // if we have a mask image and a mask planar figure, we can not do statistics // we do not know which one to use if (m_MaskImage.IsNotNull() && m_MaskPlanarFigure.IsNotNull()) { return false; } return true; } void QmitkCESTStatisticsView::Clear() { this->m_zSpectrum.clear(); this->m_zSpectrum.resize(0); this->m_ZImage = nullptr; this->m_MaskImage = nullptr; this->m_MaskPlanarFigure = nullptr; this->m_PointSet = nullptr; this->m_Controls.m_DataViewWidget->Clear(); this->InvalidateStatisticsTableView(); this->m_Controls.m_StatisticsGroupBox->setEnabled(false); } void QmitkCESTStatisticsView::OnCopyStatisticsToClipboardPushButtonClicked() { QLocale tempLocal; QLocale::setDefault(QLocale(QLocale::English, QLocale::UnitedStates)); const std::vector &statistics = this->m_CalculatorThread->GetStatisticsData(); QmitkPlotWidget::DataVector::size_type size = m_zSpectrum.size(); QString clipboard("delta_w \t Mean \t StdDev \t RMS \t Max \t Min \t N\n"); for (QmitkPlotWidget::DataVector::size_type index = 0; index < size; ++index) { // Copy statistics to clipboard ("%Ln" will use the default locale for // number formatting) clipboard = clipboard.append("%L1 \t %L2 \t %L3 \t %L4 \t %L5 \t %L6 \t %L7\n") .arg(m_zSpectrum[index], 0, 'f', 10) .arg(statistics[index]->GetMean(), 0, 'f', 10) .arg(statistics[index]->GetStd(), 0, 'f', 10) .arg(statistics[index]->GetRMS(), 0, 'f', 10) .arg(statistics[index]->GetMax(), 0, 'f', 10) .arg(statistics[index]->GetMin(), 0, 'f', 10) .arg(statistics[index]->GetN()); } QApplication::clipboard()->setText( clipboard, QClipboard::Clipboard); QLocale::setDefault(tempLocal); } void QmitkCESTStatisticsView::OnSliceChanged() { mitk::Point3D currentSelectedPosition = this->GetRenderWindowPart()->GetSelectedPosition(nullptr); unsigned int currentSelectedTimeStep = this->GetRenderWindowPart()->GetTimeNavigationController()->GetTime()->GetPos(); if (m_currentSelectedPosition != currentSelectedPosition || m_currentSelectedTimeStep != currentSelectedTimeStep) //|| m_selectedNodeTime > m_currentPositionTime) { //the current position has been changed or the selected node has been changed since the last position validation -> check position m_currentSelectedPosition = currentSelectedPosition; m_currentSelectedTimeStep = currentSelectedTimeStep; m_currentPositionTime.Modified(); m_CrosshairPointSet->Clear(); m_CrosshairPointSet->SetPoint(0, m_currentSelectedPosition); QList nodes = this->GetDataManagerSelection(); if (nodes.empty() || nodes.size() > 1) return; mitk::DataNode* node = nodes.front(); if (!node) { return; } if (dynamic_cast(node->GetData()) != nullptr) { m_Controls.labelWarning->setVisible(false); bool zSpectrumSet = SetZSpectrum(dynamic_cast( node->GetData()->GetProperty(mitk::CustomTagParser::m_OffsetsPropertyName.c_str()).GetPointer())); if (zSpectrumSet) { m_ZImage = dynamic_cast(node->GetData()); } else { return; } } else { return; } this->m_Controls.m_DataViewWidget->Clear(); AccessFixedDimensionByItk(m_ZImage, PlotPointSet, 4); } }