Medical Imaging Interaction Toolkit  2016.11.0
Medical Imaging Interaction Toolkit
mitkFiberfoxParameters.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,
6 Division of Medical and Biological Informatics.
7 All rights reserved.
8 
9 This software is distributed WITHOUT ANY WARRANTY; without
10 even the implied warranty of MERCHANTABILITY or FITNESS FOR
11 A PARTICULAR PURPOSE.
12 
13 See LICENSE.txt or http://www.mitk.org for details.
14 
15 ===================================================================*/
16 
17 #define RAPIDXML_NO_EXCEPTIONS
18 #include <boost/foreach.hpp>
19 #include <boost/lexical_cast.hpp>
20 #include <itkImageFileWriter.h>
21 #include <itkImageFileReader.h>
22 #include <mitkLog.h>
23 #include <algorithm>
24 
25 template< class ScalarType >
26 mitk::FiberfoxParameters< ScalarType >::FiberfoxParameters()
27  : m_NoiseModel(NULL)
28 {
29 
30 }
31 
32 template< class ScalarType >
33 mitk::FiberfoxParameters< ScalarType >::~FiberfoxParameters()
34 {
35  // if (m_NoiseModel!=NULL)
36  // delete m_NoiseModel;
37 }
38 
39 void mitk::SignalGenerationParameters::GenerateGradientHalfShell()
40 {
41  int NPoints = 2*m_NumGradients;
42  m_GradientDirections.clear();
43 
44  m_NumBaseline = NPoints/20;
45  if (m_NumBaseline==0)
46  m_NumBaseline=1;
47 
48  GradientType g;
49  g.Fill(0.0);
50  for (unsigned int i=0; i<m_NumBaseline; i++)
51  m_GradientDirections.push_back(g);
52 
53  if (NPoints==0)
54  return;
55 
56  vnl_vector<double> theta; theta.set_size(NPoints);
57  vnl_vector<double> phi; phi.set_size(NPoints);
58  double C = sqrt(4*M_PI);
59  phi(0) = 0.0;
60  phi(NPoints-1) = 0.0;
61  for(int i=0; i<NPoints; i++)
62  {
63  theta(i) = acos(-1.0+2.0*i/(NPoints-1.0)) - M_PI / 2.0;
64  if( i>0 && i<NPoints-1)
65  {
66  phi(i) = (phi(i-1) + C /
67  sqrt(NPoints*(1-(-1.0+2.0*i/(NPoints-1.0))*(-1.0+2.0*i/(NPoints-1.0)))));
68  }
69  }
70 
71  for(int i=0; i<NPoints; i++)
72  {
73  g[2] = sin(theta(i));
74  if (g[2]<0)
75  continue;
76  g[0] = cos(theta(i)) * cos(phi(i));
77  g[1] = cos(theta(i)) * sin(phi(i));
78  m_GradientDirections.push_back(g);
79  }
80 }
81 
82 std::vector< int > mitk::SignalGenerationParameters::GetBaselineIndices()
83 {
84  std::vector< int > result;
85  for( unsigned int i=0; i<this->m_GradientDirections.size(); i++)
86  if (m_GradientDirections.at(i).GetNorm()<0.0001)
87  result.push_back(i);
88  return result;
89 }
90 
91 unsigned int mitk::SignalGenerationParameters::GetFirstBaselineIndex()
92 {
93  for( unsigned int i=0; i<this->m_GradientDirections.size(); i++)
94  if (m_GradientDirections.at(i).GetNorm()<0.0001)
95  return i;
96  return -1;
97 }
98 
99 bool mitk::SignalGenerationParameters::IsBaselineIndex(unsigned int idx)
100 {
101  if (m_GradientDirections.size()>idx && m_GradientDirections.at(idx).GetNorm()<0.0001)
102  return true;
103  return false;
104 }
105 
106 unsigned int mitk::SignalGenerationParameters::GetNumWeightedVolumes()
107 {
108  return m_NumGradients;
109 }
110 
111 unsigned int mitk::SignalGenerationParameters::GetNumBaselineVolumes()
112 {
113  return m_NumBaseline;
114 }
115 
116 unsigned int mitk::SignalGenerationParameters::GetNumVolumes()
117 {
118  return m_GradientDirections.size();
119 }
120 
121 mitk::SignalGenerationParameters::GradientListType mitk::SignalGenerationParameters::GetGradientDirections()
122 {
123  return m_GradientDirections;
124 }
125 
126 mitk::SignalGenerationParameters::GradientType mitk::SignalGenerationParameters::GetGradientDirection(unsigned int i)
127 {
128  return m_GradientDirections.at(i);
129 }
130 
131 void mitk::SignalGenerationParameters::SetNumWeightedVolumes(int numGradients)
132 {
133  m_NumGradients = numGradients;
134  GenerateGradientHalfShell();
135 }
136 
137 void mitk::SignalGenerationParameters::SetGradienDirections(GradientListType gradientList)
138 {
139  m_GradientDirections = gradientList;
140  m_NumGradients = 0;
141  m_NumBaseline = 0;
142  for( unsigned int i=0; i<this->m_GradientDirections.size(); i++)
143  {
144  if (m_GradientDirections.at(i).GetNorm()>0.0001)
145  m_NumGradients++;
146  else
147  m_NumBaseline++;
148  }
149 }
150 
151 void mitk::SignalGenerationParameters::SetGradienDirections(mitk::DiffusionPropertyHelper::GradientDirectionsContainerType::Pointer gradientList)
152 {
153  m_NumGradients = 0;
154  m_NumBaseline = 0;
155  m_GradientDirections.clear();
156 
157  for( unsigned int i=0; i<gradientList->Size(); i++)
158  {
159  GradientType g;
160  g[0] = gradientList->at(i)[0];
161  g[1] = gradientList->at(i)[1];
162  g[2] = gradientList->at(i)[2];
163  m_GradientDirections.push_back(g);
164 
165  if (m_GradientDirections.at(i).GetNorm()>0.0001)
166  m_NumGradients++;
167  else
168  m_NumBaseline++;
169  }
170 }
171 
172 template< class ScalarType >
173 void mitk::FiberfoxParameters< ScalarType >::SaveParameters(string filename)
174 {
175  if(filename.empty())
176  return;
177  if(".ffp"!=filename.substr(filename.size()-4, 4))
178  filename += ".ffp";
179 
180  const std::string& locale = "C";
181  const std::string& currLocale = setlocale( LC_ALL, NULL );
182 
183  if ( locale.compare(currLocale)!=0 )
184  {
185  try
186  {
187  setlocale(LC_ALL, locale.c_str());
188  }
189  catch(...)
190  {
191  MITK_INFO << "Could not set locale " << locale;
192  }
193  }
194 
195  boost::property_tree::ptree parameters;
196 
197  // fiber generation parameters
198  parameters.put("fiberfox.fibers.distribution", m_FiberGen.m_Distribution);
199  parameters.put("fiberfox.fibers.variance", m_FiberGen.m_Variance);
200  parameters.put("fiberfox.fibers.density", m_FiberGen.m_Density);
201  parameters.put("fiberfox.fibers.spline.sampling", m_FiberGen.m_Sampling);
202  parameters.put("fiberfox.fibers.spline.tension", m_FiberGen.m_Tension);
203  parameters.put("fiberfox.fibers.spline.continuity", m_FiberGen.m_Continuity);
204  parameters.put("fiberfox.fibers.spline.bias", m_FiberGen.m_Bias);
205  parameters.put("fiberfox.fibers.rotation.x", m_FiberGen.m_Rotation[0]);
206  parameters.put("fiberfox.fibers.rotation.y", m_FiberGen.m_Rotation[1]);
207  parameters.put("fiberfox.fibers.rotation.z", m_FiberGen.m_Rotation[2]);
208  parameters.put("fiberfox.fibers.translation.x", m_FiberGen.m_Translation[0]);
209  parameters.put("fiberfox.fibers.translation.y", m_FiberGen.m_Translation[1]);
210  parameters.put("fiberfox.fibers.translation.z", m_FiberGen.m_Translation[2]);
211  parameters.put("fiberfox.fibers.scale.x", m_FiberGen.m_Scale[0]);
212  parameters.put("fiberfox.fibers.scale.y", m_FiberGen.m_Scale[1]);
213  parameters.put("fiberfox.fibers.scale.z", m_FiberGen.m_Scale[2]);
214 
215  // image generation parameters
216  parameters.put("fiberfox.image.basic.size.x", m_SignalGen.m_ImageRegion.GetSize(0));
217  parameters.put("fiberfox.image.basic.size.y", m_SignalGen.m_ImageRegion.GetSize(1));
218  parameters.put("fiberfox.image.basic.size.z", m_SignalGen.m_ImageRegion.GetSize(2));
219  parameters.put("fiberfox.image.basic.spacing.x", m_SignalGen.m_ImageSpacing[0]);
220  parameters.put("fiberfox.image.basic.spacing.y", m_SignalGen.m_ImageSpacing[1]);
221  parameters.put("fiberfox.image.basic.spacing.z", m_SignalGen.m_ImageSpacing[2]);
222  parameters.put("fiberfox.image.basic.origin.x", m_SignalGen.m_ImageOrigin[0]);
223  parameters.put("fiberfox.image.basic.origin.y", m_SignalGen.m_ImageOrigin[1]);
224  parameters.put("fiberfox.image.basic.origin.z", m_SignalGen.m_ImageOrigin[2]);
225  parameters.put("fiberfox.image.basic.direction.1", m_SignalGen.m_ImageDirection[0][0]);
226  parameters.put("fiberfox.image.basic.direction.2", m_SignalGen.m_ImageDirection[0][1]);
227  parameters.put("fiberfox.image.basic.direction.3", m_SignalGen.m_ImageDirection[0][2]);
228  parameters.put("fiberfox.image.basic.direction.4", m_SignalGen.m_ImageDirection[1][0]);
229  parameters.put("fiberfox.image.basic.direction.5", m_SignalGen.m_ImageDirection[1][1]);
230  parameters.put("fiberfox.image.basic.direction.6", m_SignalGen.m_ImageDirection[1][2]);
231  parameters.put("fiberfox.image.basic.direction.7", m_SignalGen.m_ImageDirection[2][0]);
232  parameters.put("fiberfox.image.basic.direction.8", m_SignalGen.m_ImageDirection[2][1]);
233  parameters.put("fiberfox.image.basic.direction.9", m_SignalGen.m_ImageDirection[2][2]);
234  parameters.put("fiberfox.image.basic.numgradients", m_SignalGen.GetNumWeightedVolumes());
235  for( unsigned int i=0; i<this->m_SignalGen.GetNumVolumes(); i++)
236  {
237  parameters.put("fiberfox.image.gradients."+boost::lexical_cast<string>(i)+".x", m_SignalGen.GetGradientDirection(i)[0]);
238  parameters.put("fiberfox.image.gradients."+boost::lexical_cast<string>(i)+".y", m_SignalGen.GetGradientDirection(i)[1]);
239  parameters.put("fiberfox.image.gradients."+boost::lexical_cast<string>(i)+".z", m_SignalGen.GetGradientDirection(i)[2]);
240  }
241 
242  parameters.put("fiberfox.image.acquisitiontype", m_SignalGen.m_AcquisitionType);
243  parameters.put("fiberfox.image.coilsensitivityprofile", m_SignalGen.m_CoilSensitivityProfile);
244  parameters.put("fiberfox.image.numberofcoils", m_SignalGen.m_NumberOfCoils);
245  parameters.put("fiberfox.image.reversephase", m_SignalGen.m_ReversePhase);
246  parameters.put("fiberfox.image.partialfourier", m_SignalGen.m_PartialFourier);
247  parameters.put("fiberfox.image.noisevariance", m_SignalGen.m_NoiseVariance);
248  parameters.put("fiberfox.image.trep", m_SignalGen.m_tRep);
249  parameters.put("fiberfox.image.signalScale", m_SignalGen.m_SignalScale);
250  parameters.put("fiberfox.image.tEcho", m_SignalGen.m_tEcho);
251  parameters.put("fiberfox.image.tLine", m_SignalGen.m_tLine);
252  parameters.put("fiberfox.image.tInhom", m_SignalGen.m_tInhom);
253  parameters.put("fiberfox.image.bvalue", m_SignalGen.m_Bvalue);
254  parameters.put("fiberfox.image.simulatekspace", m_SignalGen.m_SimulateKspaceAcquisition);
255  parameters.put("fiberfox.image.axonRadius", m_SignalGen.m_AxonRadius);
256  parameters.put("fiberfox.image.doSimulateRelaxation", m_SignalGen.m_DoSimulateRelaxation);
257  parameters.put("fiberfox.image.doDisablePartialVolume", m_SignalGen.m_DoDisablePartialVolume);
258  parameters.put("fiberfox.image.artifacts.spikesnum", m_SignalGen.m_Spikes);
259  parameters.put("fiberfox.image.artifacts.spikesscale", m_SignalGen.m_SpikeAmplitude);
260  parameters.put("fiberfox.image.artifacts.kspaceLineOffset", m_SignalGen.m_KspaceLineOffset);
261  parameters.put("fiberfox.image.artifacts.eddyStrength", m_SignalGen.m_EddyStrength);
262  parameters.put("fiberfox.image.artifacts.eddyTau", m_SignalGen.m_Tau);
263  parameters.put("fiberfox.image.artifacts.aliasingfactor", m_SignalGen.m_CroppingFactor);
264  parameters.put("fiberfox.image.artifacts.addringing", m_SignalGen.m_DoAddGibbsRinging);
265  parameters.put("fiberfox.image.artifacts.doAddMotion", m_SignalGen.m_DoAddMotion);
266  parameters.put("fiberfox.image.artifacts.randomMotion", m_SignalGen.m_DoRandomizeMotion);
267  parameters.put("fiberfox.image.artifacts.translation0", m_SignalGen.m_Translation[0]);
268  parameters.put("fiberfox.image.artifacts.translation1", m_SignalGen.m_Translation[1]);
269  parameters.put("fiberfox.image.artifacts.translation2", m_SignalGen.m_Translation[2]);
270  parameters.put("fiberfox.image.artifacts.rotation0", m_SignalGen.m_Rotation[0]);
271  parameters.put("fiberfox.image.artifacts.rotation1", m_SignalGen.m_Rotation[1]);
272  parameters.put("fiberfox.image.artifacts.rotation2", m_SignalGen.m_Rotation[2]);
273  parameters.put("fiberfox.image.artifacts.motionvolumes", m_Misc.m_MotionVolumesBox);
274  parameters.put("fiberfox.image.artifacts.addnoise", m_Misc.m_CheckAddNoiseBox);
275  parameters.put("fiberfox.image.artifacts.addghosts", m_Misc.m_CheckAddGhostsBox);
276  parameters.put("fiberfox.image.artifacts.addaliasing", m_Misc.m_CheckAddAliasingBox);
277  parameters.put("fiberfox.image.artifacts.addspikes", m_Misc.m_CheckAddSpikesBox);
278  parameters.put("fiberfox.image.artifacts.addeddycurrents", m_Misc.m_CheckAddEddyCurrentsBox);
279  parameters.put("fiberfox.image.artifacts.doAddDistortions", m_Misc.m_CheckAddDistortionsBox);
280  parameters.put("fiberfox.image.outputvolumefractions", m_Misc.m_CheckOutputVolumeFractionsBox);
281  parameters.put("fiberfox.image.showadvanced", m_Misc.m_CheckAdvancedSignalOptionsBox);
282  parameters.put("fiberfox.image.signalmodelstring", m_Misc.m_SignalModelString);
283  parameters.put("fiberfox.image.artifactmodelstring", m_Misc.m_ArtifactModelString);
284  parameters.put("fiberfox.image.outpath", m_Misc.m_OutputPath);
285  parameters.put("fiberfox.fibers.realtime", m_Misc.m_CheckRealTimeFibersBox);
286  parameters.put("fiberfox.fibers.showadvanced", m_Misc.m_CheckAdvancedFiberOptionsBox);
287  parameters.put("fiberfox.fibers.constantradius", m_Misc.m_CheckConstantRadiusBox);
288  parameters.put("fiberfox.fibers.includeFiducials", m_Misc.m_CheckIncludeFiducialsBox);
289 
290  if (m_NoiseModel!=NULL)
291  {
292  parameters.put("fiberfox.image.artifacts.noisevariance", m_NoiseModel->GetNoiseVariance());
293  if (dynamic_cast<mitk::RicianNoiseModel<ScalarType>*>(m_NoiseModel.get()))
294  parameters.put("fiberfox.image.artifacts.noisetype", "rice");
295  else if (dynamic_cast<mitk::ChiSquareNoiseModel<ScalarType>*>(m_NoiseModel.get()))
296  parameters.put("fiberfox.image.artifacts.noisetype", "chisquare");
297  }
298 
299  for (int i=0; i<m_FiberModelList.size()+m_NonFiberModelList.size(); i++)
300  {
301  mitk::DiffusionSignalModel<ScalarType>* signalModel = NULL;
302  if (i<m_FiberModelList.size())
303  {
304  signalModel = m_FiberModelList.at(i);
305  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".type", "fiber");
306  }
307  else
308  {
309  signalModel = m_NonFiberModelList.at(i-m_FiberModelList.size());
310  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".type", "non-fiber");
311  }
312 
313  if (dynamic_cast<mitk::StickModel<ScalarType>*>(signalModel))
314  {
315  mitk::StickModel<ScalarType>* model = dynamic_cast<mitk::StickModel<ScalarType>*>(signalModel);
316  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".model", "stick");
317  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".d", model->GetDiffusivity());
318  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t2", model->GetT2());
319  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t1", model->GetT1());
320  }
321  else if (dynamic_cast<mitk::TensorModel<ScalarType>*>(signalModel))
322  {
323  mitk::TensorModel<ScalarType>* model = dynamic_cast<mitk::TensorModel<ScalarType>*>(signalModel);
324  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".model", "tensor");
325  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".d1", model->GetDiffusivity1());
326  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".d2", model->GetDiffusivity2());
327  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".d3", model->GetDiffusivity3());
328  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t2", model->GetT2());
329  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t1", model->GetT1());
330  }
331  else if (dynamic_cast<mitk::RawShModel<ScalarType>*>(signalModel))
332  {
333  mitk::RawShModel<ScalarType>* model = dynamic_cast<mitk::RawShModel<ScalarType>*>(signalModel);
334  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".model", "prototype");
335  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".minFA", model->GetFaRange().first);
336  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".maxFA", model->GetFaRange().second);
337  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".minADC", model->GetAdcRange().first);
338  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".maxADC", model->GetAdcRange().second);
339  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".maxNumSamples", model->GetMaxNumKernels());
340  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".numSamples", model->GetNumberOfKernels());
341  int shOrder = model->GetShOrder();
342  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".numCoeffs", (shOrder*shOrder + shOrder + 2)/2 + shOrder);
343 
344  for (unsigned int j=0; j<model->GetNumberOfKernels(); j++)
345  {
346  vnl_vector< double > coeffs = model->GetCoefficients(j);
347  for (unsigned int k=0; k<coeffs.size(); k++)
348  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".kernels."+boost::lexical_cast<string>(j)+".coeffs."+boost::lexical_cast<string>(k), coeffs[k]);
349  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".kernels."+boost::lexical_cast<string>(j)+".B0", model->GetBaselineSignal(j));
350  }
351  }
352  else if (dynamic_cast<mitk::BallModel<ScalarType>*>(signalModel))
353  {
354  mitk::BallModel<ScalarType>* model = dynamic_cast<mitk::BallModel<ScalarType>*>(signalModel);
355  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".model", "ball");
356  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".d", model->GetDiffusivity());
357  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t2", model->GetT2());
358  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t1", model->GetT1());
359  }
360  else if (dynamic_cast<mitk::AstroStickModel<ScalarType>*>(signalModel))
361  {
362  mitk::AstroStickModel<ScalarType>* model = dynamic_cast<mitk::AstroStickModel<ScalarType>*>(signalModel);
363  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".model", "astrosticks");
364  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".d", model->GetDiffusivity());
365  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t2", model->GetT2());
366  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t1", model->GetT1());
367  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".randomize", model->GetRandomizeSticks());
368  }
369  else if (dynamic_cast<mitk::DotModel<ScalarType>*>(signalModel))
370  {
371  mitk::DotModel<ScalarType>* model = dynamic_cast<mitk::DotModel<ScalarType>*>(signalModel);
372  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".model", "dot");
373  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t2", model->GetT2());
374  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".t1", model->GetT1());
375  }
376 
377  if (signalModel!=NULL)
378  {
379  parameters.put("fiberfox.image.compartments."+boost::lexical_cast<string>(i)+".ID", signalModel->m_CompartmentId);
380 
381  if (signalModel->GetVolumeFractionImage().IsNotNull())
382  {
383  try{
384  itk::ImageFileWriter<ItkDoubleImgType>::Pointer writer = itk::ImageFileWriter<ItkDoubleImgType>::New();
385  writer->SetFileName(filename+"_VOLUME"+boost::lexical_cast<string>(signalModel->m_CompartmentId)+".nrrd");
386  writer->SetInput(signalModel->GetVolumeFractionImage());
387  writer->Update();
388  MITK_INFO << "Volume fraction image for compartment "+boost::lexical_cast<string>(signalModel->m_CompartmentId)+" saved.";
389  }
390  catch(...)
391  {
392  }
393  }
394  }
395  }
396 
397  boost::property_tree::xml_writer_settings<std::string> writerSettings(' ', 2);
398  boost::property_tree::xml_parser::write_xml(filename, parameters, std::locale(), writerSettings);
399 
400  try{
401  itk::ImageFileWriter<ItkDoubleImgType>::Pointer writer = itk::ImageFileWriter<ItkDoubleImgType>::New();
402  writer->SetFileName(filename+"_FMAP.nrrd");
403  writer->SetInput(m_SignalGen.m_FrequencyMap);
404  writer->Update();
405  }
406  catch(...)
407  {
408  MITK_INFO << "No frequency map saved.";
409  }
410  try{
411  itk::ImageFileWriter<ItkUcharImgType>::Pointer writer = itk::ImageFileWriter<ItkUcharImgType>::New();
412  writer->SetFileName(filename+"_MASK.nrrd");
413  writer->SetInput(m_SignalGen.m_MaskImage);
414  writer->Update();
415  }
416  catch(...)
417  {
418  MITK_INFO << "No mask image saved.";
419  }
420 
421  setlocale(LC_ALL, currLocale.c_str());
422 }
423 
424 
425 template< class ScalarType >
426 template< class ParameterType >
427 ParameterType mitk::FiberfoxParameters< ScalarType >::ReadVal(boost::property_tree::ptree::value_type const& v, std::string tag, ParameterType defaultValue, bool essential)
428 {
429  try
430  {
431  return v.second.get<ParameterType>(tag);
432  }
433  catch (...)
434  {
435  if (essential)
436  {
437  mitkThrow() << "Parameter file corrupted. Essential tag is missing: '" << tag << "'";
438  }
439  if (tag!="artifacts.noisetype")
440  {
441  MITK_INFO << "Tag '" << tag << "' not found. Using default value '" << defaultValue << "'.";
442  m_MissingTags += "\n- ";
443  m_MissingTags += tag;
444  }
445  return defaultValue;
446  }
447 }
448 
449 template< class ScalarType >
450 void mitk::FiberfoxParameters< ScalarType >::LoadParameters(string filename)
451 {
452  m_MissingTags = "";
453  if(filename.empty()) { return; }
454 
455  const std::string& locale = "C";
456  const std::string& currLocale = setlocale( LC_ALL, NULL );
457 
458  if ( locale.compare(currLocale)!=0 )
459  {
460  try
461  {
462  setlocale(LC_ALL, locale.c_str());
463  }
464  catch(...)
465  {
466  MITK_INFO << "Could not set locale " << locale;
467  }
468  }
469 
470 
471  boost::property_tree::ptree parameterTree;
472  boost::property_tree::xml_parser::read_xml( filename, parameterTree );
473 
474 
475  m_FiberModelList.clear();
476  m_NonFiberModelList.clear();
477  if (m_NoiseModel) { m_NoiseModel = nullptr; }
478 
479  BOOST_FOREACH( boost::property_tree::ptree::value_type const& v1, parameterTree.get_child("fiberfox") )
480  {
481  if( v1.first == "fibers" )
482  {
483  m_Misc.m_CheckRealTimeFibersBox = ReadVal<bool>(v1,"realtime", m_Misc.m_CheckRealTimeFibersBox);
484  m_Misc.m_CheckAdvancedFiberOptionsBox = ReadVal<bool>(v1,"showadvanced", m_Misc.m_CheckAdvancedFiberOptionsBox);
485  m_Misc.m_CheckConstantRadiusBox = ReadVal<bool>(v1,"constantradius", m_Misc.m_CheckConstantRadiusBox);
486  m_Misc.m_CheckIncludeFiducialsBox = ReadVal<bool>(v1,"includeFiducials", m_Misc.m_CheckIncludeFiducialsBox);
487 
488  switch (ReadVal<unsigned int>(v1,"distribution", 0))
489  {
490  case 0:
491  m_FiberGen.m_Distribution = FiberGenerationParameters::DISTRIBUTE_UNIFORM;
492  break;
493  case 1:
494  m_FiberGen.m_Distribution = FiberGenerationParameters::DISTRIBUTE_GAUSSIAN;
495  break;
496  default:
497  m_FiberGen.m_Distribution = FiberGenerationParameters::DISTRIBUTE_UNIFORM;
498  }
499  m_FiberGen.m_Variance = ReadVal<double>(v1,"variance", m_FiberGen.m_Variance);
500  m_FiberGen.m_Density = ReadVal<unsigned int>(v1,"density", m_FiberGen.m_Density);
501  m_FiberGen.m_Sampling = ReadVal<double>(v1,"spline.sampling", m_FiberGen.m_Sampling);
502  m_FiberGen.m_Tension = ReadVal<double>(v1,"spline.tension", m_FiberGen.m_Tension);
503  m_FiberGen.m_Continuity = ReadVal<double>(v1,"spline.continuity", m_FiberGen.m_Continuity);
504  m_FiberGen.m_Bias = ReadVal<double>(v1,"spline.bias", m_FiberGen.m_Bias);
505  m_FiberGen.m_Rotation[0] = ReadVal<double>(v1,"rotation.x", m_FiberGen.m_Rotation[0]);
506  m_FiberGen.m_Rotation[1] = ReadVal<double>(v1,"rotation.y", m_FiberGen.m_Rotation[1]);
507  m_FiberGen.m_Rotation[2] = ReadVal<double>(v1,"rotation.z", m_FiberGen.m_Rotation[2]);
508  m_FiberGen.m_Translation[0] = ReadVal<double>(v1,"translation.x", m_FiberGen.m_Translation[0]);
509  m_FiberGen.m_Translation[1] = ReadVal<double>(v1,"translation.y", m_FiberGen.m_Translation[1]);
510  m_FiberGen.m_Translation[2] = ReadVal<double>(v1,"translation.z", m_FiberGen.m_Translation[2]);
511  m_FiberGen.m_Scale[0] = ReadVal<double>(v1,"scale.x", m_FiberGen.m_Scale[0]);
512  m_FiberGen.m_Scale[1] = ReadVal<double>(v1,"scale.y", m_FiberGen.m_Scale[1]);
513  m_FiberGen.m_Scale[2] = ReadVal<double>(v1,"scale.z", m_FiberGen.m_Scale[2]);
514  }
515  else if ( v1.first == "image" )
516  {
517  m_Misc.m_SignalModelString = ReadVal<string>(v1,"signalmodelstring", m_Misc.m_SignalModelString);
518  m_Misc.m_ArtifactModelString = ReadVal<string>(v1,"artifactmodelstring", m_Misc.m_ArtifactModelString);
519  m_Misc.m_OutputPath = ReadVal<string>(v1,"outpath", m_Misc.m_OutputPath);
520  m_Misc.m_CheckOutputVolumeFractionsBox = ReadVal<bool>(v1,"outputvolumefractions", m_Misc.m_CheckOutputVolumeFractionsBox);
521  m_Misc.m_CheckAdvancedSignalOptionsBox = ReadVal<bool>(v1,"showadvanced", m_Misc.m_CheckAdvancedSignalOptionsBox);
522  m_Misc.m_CheckAddDistortionsBox = ReadVal<bool>(v1,"artifacts.doAddDistortions", m_Misc.m_CheckAddDistortionsBox);
523  m_Misc.m_CheckAddNoiseBox = ReadVal<bool>(v1,"artifacts.addnoise", m_Misc.m_CheckAddNoiseBox);
524  m_Misc.m_CheckAddGhostsBox = ReadVal<bool>(v1,"artifacts.addghosts", m_Misc.m_CheckAddGhostsBox);
525  m_Misc.m_CheckAddAliasingBox = ReadVal<bool>(v1,"artifacts.addaliasing", m_Misc.m_CheckAddAliasingBox);
526  m_Misc.m_CheckAddSpikesBox = ReadVal<bool>(v1,"artifacts.addspikes", m_Misc.m_CheckAddSpikesBox);
527  m_Misc.m_CheckAddEddyCurrentsBox = ReadVal<bool>(v1,"artifacts.addeddycurrents", m_Misc.m_CheckAddEddyCurrentsBox);
528 
529  m_SignalGen.m_ImageRegion.SetSize(0, ReadVal<int>(v1,"basic.size.x",m_SignalGen.m_ImageRegion.GetSize(0)));
530  m_SignalGen.m_ImageRegion.SetSize(1, ReadVal<int>(v1,"basic.size.y",m_SignalGen.m_ImageRegion.GetSize(1)));
531  m_SignalGen.m_ImageRegion.SetSize(2, ReadVal<int>(v1,"basic.size.z",m_SignalGen.m_ImageRegion.GetSize(2)));
532  m_SignalGen.m_ImageSpacing[0] = ReadVal<double>(v1,"basic.spacing.x",m_SignalGen.m_ImageSpacing[0]);
533  m_SignalGen.m_ImageSpacing[1] = ReadVal<double>(v1,"basic.spacing.y",m_SignalGen.m_ImageSpacing[1]);
534  m_SignalGen.m_ImageSpacing[2] = ReadVal<double>(v1,"basic.spacing.z",m_SignalGen.m_ImageSpacing[2]);
535  m_SignalGen.m_ImageOrigin[0] = ReadVal<double>(v1,"basic.origin.x",m_SignalGen.m_ImageOrigin[0]);
536  m_SignalGen.m_ImageOrigin[1] = ReadVal<double>(v1,"basic.origin.y",m_SignalGen.m_ImageOrigin[1]);
537  m_SignalGen.m_ImageOrigin[2] = ReadVal<double>(v1,"basic.origin.z",m_SignalGen.m_ImageOrigin[2]);
538  m_SignalGen.m_ImageDirection[0][0] = ReadVal<double>(v1,"basic.direction.1",m_SignalGen.m_ImageDirection[0][0]);
539  m_SignalGen.m_ImageDirection[0][1] = ReadVal<double>(v1,"basic.direction.2",m_SignalGen.m_ImageDirection[0][1]);
540  m_SignalGen.m_ImageDirection[0][2] = ReadVal<double>(v1,"basic.direction.3",m_SignalGen.m_ImageDirection[0][2]);
541  m_SignalGen.m_ImageDirection[1][0] = ReadVal<double>(v1,"basic.direction.4",m_SignalGen.m_ImageDirection[1][0]);
542  m_SignalGen.m_ImageDirection[1][1] = ReadVal<double>(v1,"basic.direction.5",m_SignalGen.m_ImageDirection[1][1]);
543  m_SignalGen.m_ImageDirection[1][2] = ReadVal<double>(v1,"basic.direction.6",m_SignalGen.m_ImageDirection[1][2]);
544  m_SignalGen.m_ImageDirection[2][0] = ReadVal<double>(v1,"basic.direction.7",m_SignalGen.m_ImageDirection[2][0]);
545  m_SignalGen.m_ImageDirection[2][1] = ReadVal<double>(v1,"basic.direction.8",m_SignalGen.m_ImageDirection[2][1]);
546  m_SignalGen.m_ImageDirection[2][2] = ReadVal<double>(v1,"basic.direction.9",m_SignalGen.m_ImageDirection[2][2]);
547 
548  m_SignalGen.m_AcquisitionType = (SignalGenerationParameters::AcquisitionType)
549  ReadVal<int>(v1,"acquisitiontype", m_SignalGen.m_AcquisitionType);
550  m_SignalGen.m_CoilSensitivityProfile = (SignalGenerationParameters::CoilSensitivityProfile)
551  ReadVal<int>(v1,"coilsensitivityprofile", m_SignalGen.m_CoilSensitivityProfile);
552  m_SignalGen.m_NumberOfCoils = ReadVal<int>(v1,"numberofcoils", m_SignalGen.m_NumberOfCoils);
553  m_SignalGen.m_ReversePhase = ReadVal<bool>(v1,"reversephase", m_SignalGen.m_ReversePhase);
554  m_SignalGen.m_PartialFourier = ReadVal<double>(v1,"partialfourier", m_SignalGen.m_PartialFourier);
555  m_SignalGen.m_NoiseVariance = ReadVal<double>(v1,"noisevariance", m_SignalGen.m_NoiseVariance);
556  m_SignalGen.m_tRep = ReadVal<double>(v1,"trep", m_SignalGen.m_tRep);
557  m_SignalGen.m_SignalScale = ReadVal<double>(v1,"signalScale", m_SignalGen.m_SignalScale);
558  m_SignalGen.m_tEcho = ReadVal<double>(v1,"tEcho", m_SignalGen.m_tEcho);
559  m_SignalGen.m_tLine = ReadVal<double>(v1,"tLine", m_SignalGen.m_tLine);
560  m_SignalGen.m_tInhom = ReadVal<double>(v1,"tInhom", m_SignalGen.m_tInhom);
561  m_SignalGen.m_Bvalue = ReadVal<double>(v1,"bvalue", m_SignalGen.m_Bvalue);
562  m_SignalGen.m_SimulateKspaceAcquisition = ReadVal<bool>(v1,"simulatekspace", m_SignalGen.m_SimulateKspaceAcquisition);
563 
564  m_SignalGen.m_AxonRadius = ReadVal<double>(v1,"axonRadius", m_SignalGen.m_AxonRadius);
565  m_SignalGen.m_Spikes = ReadVal<unsigned int>(v1,"artifacts.spikesnum", m_SignalGen.m_Spikes);
566  m_SignalGen.m_SpikeAmplitude = ReadVal<double>(v1,"artifacts.spikesscale", m_SignalGen.m_SpikeAmplitude);
567  m_SignalGen.m_KspaceLineOffset = ReadVal<double>(v1,"artifacts.kspaceLineOffset", m_SignalGen.m_KspaceLineOffset);
568  m_SignalGen.m_EddyStrength = ReadVal<double>(v1,"artifacts.eddyStrength", m_SignalGen.m_EddyStrength);
569  m_SignalGen.m_Tau = ReadVal<double>(v1,"artifacts.eddyTau", m_SignalGen.m_Tau);
570  m_SignalGen.m_CroppingFactor = ReadVal<double>(v1,"artifacts.aliasingfactor", m_SignalGen.m_CroppingFactor);
571  m_SignalGen.m_DoAddGibbsRinging = ReadVal<bool>(v1,"artifacts.addringing", m_SignalGen.m_DoAddGibbsRinging);
572  m_SignalGen.m_DoSimulateRelaxation = ReadVal<bool>(v1,"doSimulateRelaxation", m_SignalGen.m_DoSimulateRelaxation);
573  m_SignalGen.m_DoDisablePartialVolume = ReadVal<bool>(v1,"doDisablePartialVolume", m_SignalGen.m_DoDisablePartialVolume);
574  m_SignalGen.m_DoAddMotion = ReadVal<bool>(v1,"artifacts.doAddMotion", m_SignalGen.m_DoAddMotion);
575  m_SignalGen.m_DoRandomizeMotion = ReadVal<bool>(v1,"artifacts.randomMotion", m_SignalGen.m_DoRandomizeMotion);
576  m_SignalGen.m_Translation[0] = ReadVal<double>(v1,"artifacts.translation0", m_SignalGen.m_Translation[0]);
577  m_SignalGen.m_Translation[1] = ReadVal<double>(v1,"artifacts.translation1", m_SignalGen.m_Translation[1]);
578  m_SignalGen.m_Translation[2] = ReadVal<double>(v1,"artifacts.translation2", m_SignalGen.m_Translation[2]);
579  m_SignalGen.m_Rotation[0] = ReadVal<double>(v1,"artifacts.rotation0", m_SignalGen.m_Rotation[0]);
580  m_SignalGen.m_Rotation[1] = ReadVal<double>(v1,"artifacts.rotation1", m_SignalGen.m_Rotation[1]);
581  m_SignalGen.m_Rotation[2] = ReadVal<double>(v1,"artifacts.rotation2", m_SignalGen.m_Rotation[2]);
582 
583 
584  // m_SignalGen.SetNumWeightedVolumes(ReadVal<unsigned int>(v1,"numgradients", m_SignalGen.GetNumWeightedVolumes()));
585  SignalGenerationParameters::GradientListType gradients;
586  BOOST_FOREACH( boost::property_tree::ptree::value_type const& v2, v1.second.get_child("gradients") )
587  {
588  SignalGenerationParameters::GradientType g;
589  g[0] = ReadVal<double>(v2,"x",0);
590  g[1] = ReadVal<double>(v2,"y",0);
591  g[2] = ReadVal<double>(v2,"z",0);
592  gradients.push_back(g);
593  }
594  m_SignalGen.SetGradienDirections(gradients);
595 
596 
597  m_Misc.m_MotionVolumesBox = ReadVal<string>(v1,"artifacts.motionvolumes", m_Misc.m_MotionVolumesBox);
598  m_SignalGen.m_MotionVolumes.clear();
599 
600  if ( m_Misc.m_MotionVolumesBox == "random" )
601  {
602  for ( size_t i=0; i < m_SignalGen.GetNumVolumes(); ++i )
603  {
604  m_SignalGen.m_MotionVolumes.push_back( bool( rand()%2 ) );
605  }
606  MITK_DEBUG << "mitkFiberfoxParameters.cpp: Case m_Misc.m_MotionVolumesBox == \"random\".";
607  }
608  else if ( ! m_Misc.m_MotionVolumesBox.empty() )
609  {
610  stringstream stream( m_Misc.m_MotionVolumesBox );
611  std::vector<int> numbers;
612  int nummer = std::numeric_limits<int>::max();
613  while( stream >> nummer )
614  {
615  if( nummer < std::numeric_limits<int>::max() )
616  {
617  numbers.push_back( nummer );
618  }
619  }
620  // If a list of negative numbers is given:
621  if( *(std::min_element( numbers.begin(), numbers.end() )) < 0
622  && *(std::max_element( numbers.begin(), numbers.end() )) <= 0 ) // cave: -0 == +0
623  {
624  for ( size_t i=0; i<m_SignalGen.GetNumVolumes(); ++i )
625  {
626  m_SignalGen.m_MotionVolumes.push_back( true );
627  }
628  // set all true except those given.
629  for( auto iter = std::begin( numbers ); iter != std::end( numbers ); ++iter )
630  {
631  if ( -(*iter) < m_SignalGen.GetNumVolumes() && -(*iter) >= 0 )
632  {
633  m_SignalGen.m_MotionVolumes.at( -(*iter) ) = false;
634  }
635  }
636  MITK_DEBUG << "mitkFiberfoxParameters.cpp: Case list of negative numbers.";
637  }
638  // If a list of positive numbers is given:
639  else if( *(std::min_element( numbers.begin(), numbers.end() )) >= 0
640  && *(std::max_element( numbers.begin(), numbers.end() )) >= 0 )
641  {
642  for ( size_t i=0; i<m_SignalGen.GetNumVolumes(); ++i )
643  {
644  m_SignalGen.m_MotionVolumes.push_back( false );
645  }
646  // set all false except those given.
647  for( auto iter = std::begin( numbers ); iter != std::end( numbers ); ++iter )
648  {
649  if ( *iter < m_SignalGen.GetNumVolumes() && *iter >= 0 )
650  {
651  m_SignalGen.m_MotionVolumes.at( *iter ) = true;
652  }
653  }
654  MITK_DEBUG << "mitkFiberfoxParameters.cpp: Case list of positive numbers.";
655  }
656  else
657  {
658  MITK_WARN << "mitkFiberfoxParameters.cpp: Inconsistent list of numbers in m_MotionVolumesBox.";
659  break;
660  }
661  }
662  else
663  {
664  MITK_WARN << "mitkFiberfoxParameters.cpp: Cannot make sense of string in m_MotionVolumesBox.";
665  break;
666  }
667 
668 
669  try
670  {
671  if (ReadVal<string>(v1,"artifacts.noisetype","")=="rice")
672  {
673  m_NoiseModel = std::make_shared< mitk::RicianNoiseModel<ScalarType> >();
674  m_NoiseModel->SetNoiseVariance(ReadVal<double>(v1,"artifacts.noisevariance",m_NoiseModel->GetNoiseVariance()));
675  }
676  }
677  catch(...)
678  {
679  MITK_DEBUG << "mitkFiberfoxParameters.cpp: caught some error while trying m_NoiseModel->SetNoiseVariance()";
680  // throw;
681  }
682 
683  try
684  {
685  if (ReadVal<string>(v1,"artifacts.noisetype","")=="chisquare")
686  {
687  m_NoiseModel = std::make_shared< mitk::ChiSquareNoiseModel<ScalarType> >();
688  m_NoiseModel->SetNoiseVariance(ReadVal<double>(v1,"artifacts.noisevariance",m_NoiseModel->GetNoiseVariance()));
689  }
690  }
691  catch(...)
692  {
693  MITK_DEBUG << "mitkFiberfoxParameters.cpp: caught some error while trying m_NoiseModel->SetNoiseVariance()";
694  // throw;
695  }
696 
697 
698  BOOST_FOREACH( boost::property_tree::ptree::value_type const& v2, v1.second.get_child("compartments") )
699  {
700  mitk::DiffusionSignalModel<ScalarType>* signalModel = NULL;
701 
702  std::string model = ReadVal<std::string>(v2,"model","",true);
703  if (model=="stick")
704  {
705  mitk::StickModel<ScalarType>* model = new mitk::StickModel<ScalarType>();
706  model->SetDiffusivity(ReadVal<double>(v2,"d",model->GetDiffusivity()));
707  model->SetT2(ReadVal<double>(v2,"t2",model->GetT2()));
708  model->SetT1(ReadVal<double>(v2,"t1",model->GetT1()));
709  model->m_CompartmentId = ReadVal<unsigned int>(v2,"ID",0,true);
710  if (ReadVal<string>(v2,"type","",true)=="fiber")
711  m_FiberModelList.push_back(model);
712  else if (ReadVal<string>(v2,"type","",true)=="non-fiber")
713  m_NonFiberModelList.push_back(model);
714  signalModel = model;
715  }
716  else if (model=="tensor")
717  {
718  mitk::TensorModel<ScalarType>* model = new mitk::TensorModel<ScalarType>();
719  model->SetDiffusivity1(ReadVal<double>(v2,"d1",model->GetDiffusivity1()));
720  model->SetDiffusivity2(ReadVal<double>(v2,"d2",model->GetDiffusivity2()));
721  model->SetDiffusivity3(ReadVal<double>(v2,"d3",model->GetDiffusivity3()));
722  model->SetT2(ReadVal<double>(v2,"t2",model->GetT2()));
723  model->SetT1(ReadVal<double>(v2,"t1",model->GetT1()));
724  model->m_CompartmentId = ReadVal<unsigned int>(v2,"ID",0,true);
725  if (ReadVal<string>(v2,"type","",true)=="fiber")
726  m_FiberModelList.push_back(model);
727  else if (ReadVal<string>(v2,"type","",true)=="non-fiber")
728  m_NonFiberModelList.push_back(model);
729  signalModel = model;
730  }
731  else if (model=="ball")
732  {
733  mitk::BallModel<ScalarType>* model = new mitk::BallModel<ScalarType>();
734  model->SetDiffusivity(ReadVal<double>(v2,"d",model->GetDiffusivity()));
735  model->SetT2(ReadVal<double>(v2,"t2",model->GetT2()));
736  model->SetT1(ReadVal<double>(v2,"t1",model->GetT1()));
737  model->m_CompartmentId = ReadVal<unsigned int>(v2,"ID",0,true);
738  if (ReadVal<string>(v2,"type","",true)=="fiber")
739  m_FiberModelList.push_back(model);
740  else if (ReadVal<string>(v2,"type","",true)=="non-fiber")
741  m_NonFiberModelList.push_back(model);
742  signalModel = model;
743  }
744  else if (model=="astrosticks")
745  {
746  mitk::AstroStickModel<ScalarType>* model = new AstroStickModel<ScalarType>();
747  model->SetDiffusivity(ReadVal<double>(v2,"d",model->GetDiffusivity()));
748  model->SetT2(ReadVal<double>(v2,"t2",model->GetT2()));
749  model->SetT1(ReadVal<double>(v2,"t1",model->GetT1()));
750  model->SetRandomizeSticks(ReadVal<bool>(v2,"randomize",model->GetRandomizeSticks()));
751  model->m_CompartmentId = ReadVal<unsigned int>(v2,"ID",0,true);
752  if (ReadVal<string>(v2,"type","",true)=="fiber")
753  m_FiberModelList.push_back(model);
754  else if (ReadVal<string>(v2,"type","",true)=="non-fiber")
755  m_NonFiberModelList.push_back(model);
756  signalModel = model;
757  }
758  else if (model=="dot")
759  {
760  mitk::DotModel<ScalarType>* model = new mitk::DotModel<ScalarType>();
761  model->SetT2(ReadVal<double>(v2,"t2",model->GetT2()));
762  model->SetT1(ReadVal<double>(v2,"t1",model->GetT1()));
763  model->m_CompartmentId = ReadVal<unsigned int>(v2,"ID",0,true);
764  if (ReadVal<string>(v2,"type","",true)=="fiber")
765  m_FiberModelList.push_back(model);
766  else if (ReadVal<string>(v2,"type","",true)=="non-fiber")
767  m_NonFiberModelList.push_back(model);
768  signalModel = model;
769  }
770  else if (model=="prototype")
771  {
772  mitk::RawShModel<ScalarType>* model = new mitk::RawShModel<ScalarType>();
773  model->SetMaxNumKernels(ReadVal<unsigned int>(v2,"maxNumSamples",model->GetMaxNumKernels()));
774  model->SetFaRange(ReadVal<double>(v2,"minFA",model->GetFaRange().first), ReadVal<double>(v2,"maxFA",model->GetFaRange().second));
775  model->SetAdcRange(ReadVal<double>(v2,"minADC",model->GetAdcRange().first), ReadVal<double>(v2,"maxADC",model->GetAdcRange().second));
776  model->m_CompartmentId = ReadVal<unsigned int>(v2,"ID",0,true);
777  unsigned int numCoeffs = ReadVal<unsigned int>(v2,"numCoeffs",0,true);
778  unsigned int numSamples = ReadVal<unsigned int>(v2,"numSamples",0,true);
779 
780  for (unsigned int j=0; j<numSamples; j++)
781  {
782  vnl_vector< double > coeffs(numCoeffs);
783  for (unsigned int k=0; k<numCoeffs; k++)
784  {
785  coeffs[k] = ReadVal<double>(v2,"kernels."+boost::lexical_cast<string>(j)+".coeffs."+boost::lexical_cast<string>(k),0,true);
786  }
787  model->SetShCoefficients( coeffs, ReadVal<double>(v2,"kernels."+boost::lexical_cast<string>(j)+".B0",0,true) );
788  }
789 
790  if (ReadVal<string>(v2,"type","",true)=="fiber")
791  { m_FiberModelList.push_back(model); }
792  else if (ReadVal<string>(v2,"type","",true)=="non-fiber")
793  { m_NonFiberModelList.push_back(model); }
794  // else ?
795  signalModel = model;
796  }
797 
798  if (signalModel!=NULL)
799  {
800  signalModel->SetGradientList(gradients);
801  try
802  {
803  itk::ImageFileReader<ItkDoubleImgType>::Pointer reader = itk::ImageFileReader<ItkDoubleImgType>::New();
804  reader->SetFileName(filename+"_VOLUME"+ReadVal<string>(v2,"ID","")+".nrrd");
805  reader->Update();
806  signalModel->SetVolumeFractionImage(reader->GetOutput());
807  MITK_INFO << "Volume fraction image loaded for compartment " << signalModel->m_CompartmentId;
808  }
809  catch(...)
810  {
811  MITK_INFO << "No volume fraction image found for compartment " << signalModel->m_CompartmentId;
812  }
813  }
814  }
815  }
816  else
817  {
818 
819  }
820  }
821 
822  try
823  {
824  itk::ImageFileReader<ItkDoubleImgType>::Pointer reader = itk::ImageFileReader<ItkDoubleImgType>::New();
825  reader->SetFileName(filename+"_FMAP.nrrd");
826  reader->Update();
827  m_SignalGen.m_FrequencyMap = reader->GetOutput();
828  MITK_INFO << "Frequency map loaded.";
829  }
830  catch(...)
831  {
832  MITK_INFO << "No frequency map found.";
833  }
834 
835  try
836  {
837  itk::ImageFileReader<ItkUcharImgType>::Pointer reader = itk::ImageFileReader<ItkUcharImgType>::New();
838  reader->SetFileName(filename+"_MASK.nrrd");
839  reader->Update();
840  m_SignalGen.m_MaskImage = reader->GetOutput();
841  MITK_INFO << "Mask image loaded.";
842  }
843  catch(...)
844  {
845  MITK_INFO << "No mask image found.";
846  }
847 
848  setlocale(LC_ALL, currLocale.c_str());
849 }
850 
851 template< class ScalarType >
852 void mitk::FiberfoxParameters< ScalarType >::PrintSelf()
853 {
854  MITK_INFO << "Not implemented :(";
855 }
mitk::RawShModel::SetShCoefficients
bool SetShCoefficients(vnl_vector< double > shCoefficients, double b0)
Definition: mitkRawShModel.cpp:249
mitk::TensorModel::GetDiffusivity2
double GetDiffusivity2()
Definition: mitkTensorModel.h:66
mitk::FiberfoxParameters::ReadVal
ParameterType ReadVal(boost::property_tree::ptree::value_type const &v, std::string tag, ParameterType defaultValue, bool essential=false)
Definition: mitkFiberfoxParameters.cpp:427
mitk::Pointer
itk::SmartPointer< Self > Pointer
Definition: mitkRenderingManager.h:389
mitk::TensorModel::GetDiffusivity1
double GetDiffusivity1()
Definition: mitkTensorModel.h:65
mitk::DiffusionSignalModel::m_CompartmentId
unsigned int m_CompartmentId
GUI flag. Which compartment is this model assigned to?
Definition: mitkDiffusionSignalModel.h:81
mitk::RawShModel::SetFaRange
void SetFaRange(double min, double max)
Definition: mitkRawShModel.h:73
mitk::TensorModel
Generates diffusion measurement employing a second rank tensor model: e^(-bg^TDg) ...
Definition: mitkTensorModel.h:31
MITK_INFO
#define MITK_INFO
Definition: mitkLogMacros.h:22
mitk::BallModel::SetDiffusivity
void SetDiffusivity(double D)
Definition: mitkBallModel.h:58
mitk::SignalGenerationParameters::SetGradienDirections
void SetGradienDirections(GradientListType gradientList)
Definition: mitkFiberfoxParameters.cpp:137
mitk::FiberfoxParameters::SaveParameters
void SaveParameters(string filename)
Save image generation parameters to .ffp file.
Definition: mitkFiberfoxParameters.cpp:173
MITK_DEBUG
#define MITK_DEBUG
Definition: mitkLogMacros.h:26
mitk::DotModel
Generates constant direction independent signal.
Definition: mitkDotModel.h:30
mitk::SignalGenerationParameters::GradientType
itk::Vector< double, 3 > GradientType
Definition: mitkFiberfoxParameters.h:51
mitk::SignalGenerationParameters::GradientListType
std::vector< GradientType > GradientListType
Definition: mitkFiberfoxParameters.h:52
mitk::RawShModel::GetNumberOfKernels
unsigned int GetNumberOfKernels()
Definition: mitkRawShModel.cpp:67
mitk::BallModel
Generates direction independent diffusion measurement employing a scalar diffusion constant d: e^(-bd...
Definition: mitkBallModel.h:30
mitk::RawShModel::GetMaxNumKernels
unsigned int GetMaxNumKernels()
Definition: mitkRawShModel.h:79
mitk::RawShModel::SetAdcRange
void SetAdcRange(double min, double max)
Definition: mitkRawShModel.h:74
mitk::DiffusionSignalModel
Abstract class for diffusion signal models.
Definition: mitkDiffusionSignalModel.h:35
mitk::TensorModel::SetDiffusivity3
void SetDiffusivity3(double d3)
Definition: mitkTensorModel.h:64
mitk::RawShModel::GetBaselineSignal
double GetBaselineSignal(int index)
Definition: mitkRawShModel.h:88
mitk::SignalGenerationParameters::GetFirstBaselineIndex
unsigned int GetFirstBaselineIndex()
Returns index of first non-diffusion-weighted image volume.
Definition: mitkFiberfoxParameters.cpp:91
mitk::SignalGenerationParameters::GenerateGradientHalfShell
void GenerateGradientHalfShell()
Generates half shell of gradient directions (with m_NumGradients non-zero directions) ...
Definition: mitkFiberfoxParameters.cpp:39
mitk::TensorModel::GetDiffusivity3
double GetDiffusivity3()
Definition: mitkTensorModel.h:67
mitk::SignalGenerationParameters::GetBaselineIndices
std::vector< int > GetBaselineIndices()
Returns list of nun-diffusion-weighted image volume indices.
Definition: mitkFiberfoxParameters.cpp:82
mitk::SignalGenerationParameters::GetNumBaselineVolumes
unsigned int GetNumBaselineVolumes()
Get number of non-diffusion-weighted image volumes.
Definition: mitkFiberfoxParameters.cpp:111
MITK_WARN
#define MITK_WARN
Definition: mitkLogMacros.h:23
mitk::RawShModel
The spherical harmonic representation of a prototype diffusion weighted MR signal is used to obtain t...
Definition: mitkRawShModel.h:32
mitk::RawShModel::GetShOrder
unsigned int GetShOrder()
Definition: mitkRawShModel.h:85
mitk::RawShModel::GetFaRange
std::pair< double, double > GetFaRange()
Definition: mitkRawShModel.h:77
mitk::DiffusionSignalModel::SetGradientList
virtual void SetGradientList(GradientListType gradientList)=0
mitk::StickModel::SetDiffusivity
void SetDiffusivity(double diffusivity)
Scalar diffusion constant.
Definition: mitkStickModel.h:59
mitk::RawShModel::GetAdcRange
std::pair< double, double > GetAdcRange()
Definition: mitkRawShModel.h:78
mitk::TensorModel::SetDiffusivity1
void SetDiffusivity1(double d1)
Definition: mitkTensorModel.h:62
filename
static const std::string filename
Definition: mitkTinyXMLTest.cpp:28
mitkThrow
#define mitkThrow()
Definition: mitkExceptionMacro.h:31
mitk::SignalGenerationParameters::GetGradientDirections
GradientListType GetGradientDirections()
Return gradient direction container.
Definition: mitkFiberfoxParameters.cpp:121
mitk::FiberfoxParameters::PrintSelf
void PrintSelf()
Print parameters to stdout.
Definition: mitkFiberfoxParameters.cpp:852
mitk::AstroStickModel::SetDiffusivity
void SetDiffusivity(double diffusivity)
Scalar diffusion constant.
Definition: mitkAstroStickModel.h:71
mitk::ChiSquareNoiseModel
Implementation of noise following a chi-squared distribution.
Definition: mitkChiSquareNoiseModel.h:31
max
static T max(T x, T y)
Definition: svm.cpp:70
mitk::SignalGenerationParameters::SetNumWeightedVolumes
void SetNumWeightedVolumes(int numGradients)
Automaticall calls GenerateGradientHalfShell() afterwards.
Definition: mitkFiberfoxParameters.cpp:131
mitk::TensorModel::SetDiffusivity2
void SetDiffusivity2(double d2)
Definition: mitkTensorModel.h:63
mitk::StickModel
Generates the diffusion signal using an idealised cylinder with zero radius: e^(-bd(ng)²) ...
Definition: mitkStickModel.h:30
mitk::SignalGenerationParameters::GetNumWeightedVolumes
unsigned int GetNumWeightedVolumes()
Get number of diffusion-weighted image volumes.
Definition: mitkFiberfoxParameters.cpp:106
mitk::SignalGenerationParameters::IsBaselineIndex
bool IsBaselineIndex(unsigned int idx)
Checks if image volume with given index is non-diffusion-weighted volume or not.
Definition: mitkFiberfoxParameters.cpp:99
mitk::RicianNoiseModel
Implementation of noise following a rician distribution.
Definition: mitkRicianNoiseModel.h:31
mitk::RawShModel::SetMaxNumKernels
void SetMaxNumKernels(unsigned int max)
Definition: mitkRawShModel.h:75
mitk::DiffusionSignalModel::GetVolumeFractionImage
ItkDoubleImgType::Pointer GetVolumeFractionImage()
Definition: mitkDiffusionSignalModel.h:69
mitk::FiberfoxParameters::LoadParameters
void LoadParameters(string filename)
Load image generation parameters from .ffp file.
Definition: mitkFiberfoxParameters.cpp:450
mitk::AstroStickModel
Generates the diffusion signal using a collection of idealised cylinder with zero radius: e^(-bd(ng)²...
Definition: mitkAstroStickModel.h:31
mitk::RawShModel::GetCoefficients
vnl_vector< double > GetCoefficients(int listIndex)
Definition: mitkRawShModel.h:89
mitk::BallModel::GetDiffusivity
double GetDiffusivity()
Definition: mitkBallModel.h:59
mitk::SignalGenerationParameters::GetNumVolumes
unsigned int GetNumVolumes()
Get number of baseline and diffusion-weighted image volumes.
Definition: mitkFiberfoxParameters.cpp:116
mitk::StickModel::GetDiffusivity
double GetDiffusivity()
Definition: mitkStickModel.h:60
mitk::AstroStickModel::SetRandomizeSticks
void SetRandomizeSticks(bool randomize=true)
Random stick configuration in each voxel.
Definition: mitkAstroStickModel.h:65
mitk::SignalGenerationParameters::GetGradientDirection
GradientType GetGradientDirection(unsigned int i)
Definition: mitkFiberfoxParameters.cpp:126
mitk::DiffusionSignalModel::SetVolumeFractionImage
void SetVolumeFractionImage(ItkDoubleImgType::Pointer img)
Definition: mitkDiffusionSignalModel.h:68
mitk::New
static itkEventMacro(BoundingShapeInteractionEvent, itk::AnyEvent) class MITKBOUNDINGSHAPE_EXPORT BoundingShapeInteractor Pointer New()
Basic interaction methods for mitk::GeometryData.