Medical Imaging Interaction Toolkit  2018.4.99-12ad79a3
Medical Imaging Interaction Toolkit
mitkModelFitPlotDataHelper.cpp
Go to the documentation of this file.
1 /*============================================================================
2 
3 The Medical Imaging Interaction Toolkit (MITK)
4 
5 Copyright (c) German Cancer Research Center (DKFZ)
6 All rights reserved.
7 
8 Use of this source code is governed by a 3-clause BSD license that can be
9 found in the LICENSE file.
10 
11 ============================================================================*/
12 
14 
15 #include "mitkExceptionMacro.h"
16 #include "mitkImage.h"
18 #include "mitkModelGenerator.h"
19 
20 #include "mitkFormulaParser.h"
21 
23 {
24  return "Sample";
25 };
26 
28 {
29  return "Signal";
30 }
31 
33 {
34  return "INTERP_Signal";
35 };
36 
37 void
40 {
41  if (this->m_Values != _arg)
42  {
43  this->m_Values = _arg;
44  this->Modified();
45  }
46 }
47 
48 void
51 {
52  if (this->m_Values != _arg)
53  {
54  this->m_Values = std::move(_arg);
55  this->Modified();
56  }
57 }
58 
60 {
61  this->m_Values = rhs.m_Values;
62  this->SetTimeStamp(rhs.GetTimeStamp());
63  return *this;
64 };
65 
67 {
68  this->m_Values = std::move(rhs.m_Values);
69  this->SetTimeStamp(rhs.GetTimeStamp());
70 
71  return *this;
72 };
73 
75 {
76  this->m_Values.clear();
77  this->Modified();
78 }
79 
81 {
82 }
83 
84 const mitk::PlotDataCurve* GetPlotCurve(const mitk::PlotDataCurveCollection* collection, const std::string& key)
85 {
86  if (collection)
87  {
88  auto iter = collection->find(key);
89  if (iter != collection->end())
90  {
91  return iter->second.GetPointer();
92  }
93  }
94  return nullptr;
95 };
96 
98 {
99  if (coll)
100  {
102  }
103  return nullptr;
104 };
105 
107 {
108  if (coll)
109  {
111  }
112  return nullptr;
113 };
114 
116 {
117  if (coll)
118  {
120  }
121  return nullptr;
122 };
123 
124 std::string mitk::ModelFitPlotData::GetPositionalCollectionName(const PositionalCollectionMap::value_type& mapValue)
125 {
126  std::ostringstream nameStrm;
127  nameStrm.imbue(std::locale("C"));
128  nameStrm << "Pos " << mapValue.first << std::endl << std::setprecision(3) << "(" << mapValue.second.first[0] << "|" << mapValue.second.first[1] << "|" << mapValue.second.first[2] << ")";
129  return nameStrm.str();
130 };
131 
132 
134 {
135  auto predicate = [point](const PositionalCollectionMap::value_type& value) {return value.second.first == point; };
136 
137  auto iter = std::find_if(std::begin(this->positionalPlots), std::end(this->positionalPlots), predicate);
138  if (iter != positionalPlots.end())
139  {
140  return iter->second.second.GetPointer();
141  }
142  return nullptr;
143 };
144 
146 {
147  auto iter = this->positionalPlots.find(id);
148  if (iter != positionalPlots.end())
149  {
150  return iter->second.second.GetPointer();
151  }
152  return nullptr;
153 };
154 
155 
156 mitk::PlotDataValues::value_type mitk::ModelFitPlotData::GetXMinMax() const
157 {
158  double max = itk::NumericTraits<double>::NonpositiveMin();
160 
161  //currently we assume that within a model fit, plot data does not exceed
162  //the sample/signale on the x axis.
163  auto sample = this->GetSamplePlot(this->currentPositionPlots);
164  if (sample)
165  {
166  CheckXMinMaxFromPlotDataValues(sample->GetValues(), min, max);
167  }
168  for (const auto& posCollection : this->positionalPlots)
169  {
170  auto sample = this->GetSamplePlot(posCollection.second.second);
171  if (sample)
172  {
173  CheckXMinMaxFromPlotDataValues(sample->GetValues(), min, max);
174  }
175  }
176 
177  return std::make_pair(min, max);
178 };
179 
180 mitk::PlotDataValues::value_type mitk::ModelFitPlotData::GetYMinMax() const
181 {
182  double max = itk::NumericTraits<double>::NonpositiveMin();
184 
185  for (const auto& plot : *(this->currentPositionPlots.GetPointer()))
186  {
187  CheckYMinMaxFromPlotDataValues(plot.second->GetValues(), min, max);
188  }
189 
190  for (const auto& posCollection : this->positionalPlots)
191  {
192  for (const auto& plot : *(posCollection.second.second))
193  {
194  CheckYMinMaxFromPlotDataValues(plot.second->GetValues(), min, max);
195  }
196  }
197 
198  for (const auto& plot : *(this->staticPlots))
199  {
200  CheckYMinMaxFromPlotDataValues(plot.second->GetValues(), min, max);
201  }
202 
203  return std::make_pair(min, max);
204 };
205 
207 {
208  this->currentPositionPlots = PlotDataCurveCollection::New();
209  this->staticPlots = PlotDataCurveCollection::New();
210 };
211 
212 void mitk::CheckYMinMaxFromPlotDataValues(const PlotDataValues& data, double& min, double& max)
213 {
214  for (const auto & pos : data)
215  {
216  if (max < pos.second)
217  {
218  max = pos.second;
219  }
220 
221  if (min > pos.second)
222  {
223  min = pos.second;
224  }
225  }
226 }
227 
228 void mitk::CheckXMinMaxFromPlotDataValues(const PlotDataValues& data, double& min, double& max)
229 {
230  for (const auto & pos : data)
231  {
232  if (max < pos.first)
233  {
234  max = pos.first;
235  }
236 
237  if (min > pos.first)
238  {
239  min = pos.first;
240  }
241  }
242 }
243 
245 mitk::PlotDataCurve::Pointer
247 {
248  if (!fitInfo)
249  {
250  mitkThrow() << "Cannot calc model curve from function for given fit. Passed ModelFitInfo instance is nullptr.";
251  }
252 
253  mitk::Image::Pointer inputImage = fitInfo->inputImage;
254  assert(inputImage.IsNotNull());
255 
256  mitk::PlotDataCurve::Pointer result = mitk::PlotDataCurve::New();
258  values.reserve(timeGrid.size());
259 
260  // Calculate index
261  ::itk::Index<3> index;
262  mitk::BaseGeometry::Pointer geometry = inputImage->GetTimeGeometry()->GetGeometryForTimeStep(0);
263 
264  geometry->WorldToIndex(position, index);
265 
267 
268  mitk::FormulaParser parser(&parameterMap);
269 
270  for (unsigned int t = 0; t < timeGrid.size(); ++t)
271  {
272  // Set up static parameters
273  for (const auto& var : fitInfo->staticParamMap)
274  {
275  const auto& list = var.second;
276 
277  if (list.size() == 1)
278  {
279  parameterMap[var.first] = list.front();
280  }
281  else
282  {
283  parameterMap[var.first] = list.at(t);
284  }
285  }
286 
287  // Calculate curve data
288  double x = timeGrid[t];
289  parameterMap[fitInfo->x] = x;
290 
291  double y = parser.parse(fitInfo->function);
292  values.emplace_back(std::make_pair(x, y));
293  }
294 
295  result->SetValues(std::move(values));
296  return result;
297 }
298 
300 mitk::PlotDataCurve::Pointer
301 CalcSignalFromModel(const mitk::Point3D& position, const mitk::modelFit::ModelFitInfo* fitInfo, const mitk::ModelParameterizerBase* parameterizer = nullptr)
302 {
303  assert(fitInfo);
304 
305  if (!parameterizer)
306  {
307  parameterizer = mitk::ModelGenerator::GenerateModelParameterizer(*fitInfo);
308  }
309 
310  mitk::Image::Pointer inputImage = fitInfo->inputImage;
311  assert(inputImage.IsNotNull());
312 
313  // Calculate index
314  ::itk::Index<3> index;
315  mitk::BaseGeometry::Pointer geometry = inputImage->GetTimeGeometry()->GetGeometryForTimeStep(0);
316 
317  geometry->WorldToIndex(position, index);
318 
319  //model generation
320  mitk::ModelBase::Pointer model = parameterizer->GenerateParameterizedModel(index);
321 
323 
325 
326  mitk::ModelBase::ModelResultType curveDataY = model->GetSignal(paramArray);
327  mitk::PlotDataCurve::Pointer result = mitk::PlotDataCurve::New();
328 
329  mitk::ModelBase::TimeGridType timeGrid = model->GetTimeGrid();
331  values.reserve(timeGrid.size());
332 
333  for (unsigned int t = 0; t < timeGrid.size(); ++t)
334  {
335  double x = timeGrid[t];
336  double y = curveDataY[t];
337  values.emplace_back(std::make_pair(x, y));
338  }
339 
340  result->SetValues(std::move(values));
341  return result;
342 }
343 
344 mitk::PlotDataCurve::Pointer
346 {
347  if (!fitInfo)
348  {
349  mitkThrow() << "Cannot calc model curve from function for given fit. Passed ModelFitInfo instance is nullptr.";
350  }
351 
352  mitk::PlotDataCurve::Pointer result;
353 
354  if (!parameterizer)
355  {
356  parameterizer = ModelGenerator::GenerateModelParameterizer(*fitInfo);
357  }
358 
359  if (parameterizer)
360  {
361  // Use model instead of formula parser
362  parameterizer->SetDefaultTimeGrid(timeGrid);
363  result = CalcSignalFromModel(position, fitInfo, parameterizer);
364  }
365  else
366  {
367  // Use formula parser to parse function string
368  try
369  {
370  result = CalcSignalFromFunction(position, fitInfo, timeGrid);
371  }
372  catch (const mitk::FormulaParserException& e)
373  {
374  MITK_ERROR << "Error while parsing modelfit function: " << e.what();
375  }
376  }
377 
378  return result;
379 }
380 
381 mitk::PlotDataCurveCollection::Pointer
383 {
384  if (!fitInfo)
385  {
386  mitkThrow() << "Cannot calc model curve from function for given fit. Passed ModelFitInfo instance is nullptr.";
387  }
388 
389  mitk::PlotDataCurveCollection::Pointer result = mitk::PlotDataCurveCollection::New();
390 
391  for (const auto& additionalInput : fitInfo->inputData.GetLookupTable())
392  {
393  if (additionalInput.second.size() != timeGrid.size())
394  {
395  MITK_ERROR <<
396  "Error while refreshing input data for visualization. Size of data and input image time grid differ. Invalid data name: "
397  << additionalInput.first;
398  }
399  else
400  {
401  mitk::PlotDataCurve::Pointer pointData = mitk::PlotDataCurve::New();;
403  values.reserve(timeGrid.size());
404 
405  for (unsigned int t = 0; t < timeGrid.size(); ++t)
406  {
407  const double x = timeGrid[t];
408  const double y = additionalInput.second[t];
409  values.emplace_back(std::make_pair(x, y));
410  }
411  pointData->SetValues(std::move(values));
412  result->CastToSTLContainer().emplace(additionalInput.first, std::move(pointData));
413  }
414  }
415 
416  return result;
417 }
418 
419 mitk::PlotDataCurve::Pointer
421 {
422  if (!image)
423  {
424  mitkThrow() << "Cannot generate sample plot data. Passed image instance is nullptr.";
425  }
426 
427  mitk::PlotDataCurve::Pointer result = mitk::PlotDataCurve::New();
429  values.reserve(timeGrid.size());
430 
431  for (unsigned int t = 0; t < timeGrid.size(); ++t)
432  {
433  const double x = timeGrid[t];
434  const double y = ReadVoxel(image, position, t);
435  values.emplace_back(std::make_pair(x, y));
436  }
437 
438  result->SetValues(std::move(values));
439  return result;
440 }
mitk::PlotDataCurve::Pointer CalcSignalFromFunction(const mitk::Point3D &position, const mitk::modelFit::ModelFitInfo *fitInfo, const mitk::ModelBase::TimeGridType &timeGrid)
MITKMODELFIT_EXPORT const std::string MODEL_FIT_PLOT_SAMPLE_NAME()
MITKMODELFIT_EXPORT ModelTraitsInterface::ParameterValueType ReadVoxel(const mitk::Image *image, const mitk::Point3D &position, unsigned int timestep=0, bool noThrow=true)
ModelTraitsInterface::ParametersType ParametersType
Definition: mitkModelBase.h:59
MITKMODELFIT_EXPORT ModelTraitsInterface::ParametersType ConvertParameterMapToParameterVector(const ParameterValueMapType &valueMap, const ModelTraitsInterface *pTraitInterface)
virtual void SetDefaultTimeGrid(TimeGridType _arg)
static Pointer New()
This class offers the functionality to evaluate simple mathematical formula strings (e...
#define MITK_ERROR
Definition: mitkLogMacros.h:20
itk::Array< double > TimeGridType
Definition: mitkModelBase.h:62
mitk::Image::Pointer inputImage
static std::string GetPositionalCollectionName(const PositionalCollectionMap::value_type &mapValue)
StaticParameterMap staticParamMap
Exception class for all exceptions that are generated in the FormulaParser module.
static ModelParameterizerBase::Pointer GenerateModelParameterizer(const modelFit::ModelFitInfo &fit)
DataType::PointIdentifier PointIdentifier
Definition: mitkPointSet.h:131
std::map< ModelTraitsInterface::ParameterNameType, double > ParameterValueMapType
ValueType parse(const std::string &input)
Evaluates the input string and returns the resulting value.
const mitk::PlotDataCurve * GetPlotCurve(const mitk::PlotDataCurveCollection *collection, const std::string &key)
MITKMODELFIT_EXPORT PlotDataCurve::Pointer GenerateImageSamplePlotData(const mitk::Point3D &position, const mitk::Image *image, const mitk::ModelBase::TimeGridType &timeGrid)
MITKMODELFIT_EXPORT const std::string MODEL_FIT_PLOT_INTERPOLATED_SIGNAL_NAME()
MITKMODELFIT_EXPORT PlotDataCurveCollection::Pointer GenerateAdditionalModelFitPlotData(const mitk::Point3D &position, const mitk::modelFit::ModelFitInfo *fitInfo, const mitk::ModelBase::TimeGridType &timeGrid)
std::vector< std::pair< double, double > > PlotDataValues
MITKMODELFIT_EXPORT PlotDataCurve::Pointer GenerateModelSignalPlotData(const mitk::Point3D &position, const mitk::modelFit::ModelFitInfo *fitInfo, const mitk::ModelBase::TimeGridType &timeGrid, mitk::ModelParameterizerBase *parameterizer=nullptr)
#define mitkThrow()
Image class for storing images.
Definition: mitkImage.h:72
ModelTraitsInterface::ModelResultType ModelResultType
Definition: mitkModelBase.h:55
static T max(T x, T y)
Definition: svm.cpp:56
void CheckXMinMaxFromPlotDataValues(const PlotDataValues &data, double &min, double &max)
mitk::Image::Pointer image
Data class that stores all information about a modelfit that is relevant to the visualization and sto...
const LookupTableType & GetLookupTable() const
Returns the map of lists.
const PlotDataCurveCollection * GetPositionalPlot(const mitk::Point3D &point) const
void CheckYMinMaxFromPlotDataValues(const PlotDataValues &data, double &min, double &max)
static T min(T x, T y)
Definition: svm.cpp:53
itk::MapContainer< std::string, PlotDataCurve::Pointer > PlotDataCurveCollection
MITKMODELFIT_EXPORT const std::string MODEL_FIT_PLOT_SIGNAL_NAME()
PlotDataValues::value_type GetYMinMax() const
MITKMODELFIT_EXPORT ParameterValueMapType ExtractParameterValueMapFromModelFit(const mitk::modelFit::ModelFitInfo *fitInfo, const mitk::Point3D &position)
ScalarListLookupTable inputData
PlotDataValues::value_type GetXMinMax() const
static const PlotDataCurve * GetInterpolatedSignalPlot(const PlotDataCurveCollection *coll)
static const PlotDataCurve * GetSamplePlot(const PlotDataCurveCollection *coll)
PlotDataCurve & operator=(const PlotDataCurve &rhs)
mitk::PlotDataCurve::Pointer CalcSignalFromModel(const mitk::Point3D &position, const mitk::modelFit::ModelFitInfo *fitInfo, const mitk::ModelParameterizerBase *parameterizer=nullptr)
static const PlotDataCurve * GetSignalPlot(const PlotDataCurveCollection *coll)
virtual void SetValues(const ValuesType &_arg)