+++ /dev/null
-*.swp
-build
typedef unsigned char uchar;
typedef unsigned short ushort;
typedef unsigned int uint;
+
+#define CLITK_TRY_CATCH_EXIT(func) \
+ try { \
+ func; \
+ } \
+ catch (const itk::ExceptionObject& e) { \
+ e.Print(std::cout); \
+ exit(-1);\
+ } \
+ catch (const std::exception& e) { \
+ std::cout << e.what() << std::endl; \
+ exit(-2);\
+ } \
+ catch (...) { \
+ std::cout << "Unknown excpetion" << std::endl; \
+ exit(-3); \
+ }
+
//--------------------------------------------------------------------
// when everything goes wrong
//--------------------------------------------------------------------
void PrintMemoryUsed();
+ //--------------------------------------------------------------------
+ // Convert a map to a vector
+ template <typename M, typename V>
+ void MapToVecFirst(const M & m, V & v);
+ template <typename M, typename V>
+ void MapToVecSecond(const M & m, V & v);
+
#include "clitkCommon.txx"
} // end namespace
}
//--------------------------------------------------------------------
+
//--------------------------------------------------------------------
template<class ImageType>
void CloneImage(const typename ImageType::Pointer & input, typename ImageType::Pointer & output)
}
//--------------------------------------------------------------------
+
+//--------------------------------------------------------------------
+// http://stackoverflow.com/questions/771453/copy-map-values-to-vector-in-stl
+template <typename M, typename V>
+void MapToVecFirst(const M & m, V & v) {
+ for( typename M::const_iterator it = m.begin(); it != m.end(); ++it ) {
+ v.push_back( it->first );
+ }
+}
+template <typename M, typename V>
+void MapToVecSecond(const M & m, V & v) {
+ for( typename M::const_iterator it = m.begin(); it != m.end(); ++it ) {
+ v.push_back( it->second );
+ }
+}
+//--------------------------------------------------------------------
+
#endif /* end #define CLITKCOMMON_TXX */
//--------------------------------------------------------------------
/*
- Convenient class to manage options from GGO (gengetopt) to filter
+ Convenient class to manage frequent options
*/
//--------------------------------------------------------------------
class FilterBase
- BSD See included LICENSE.txt file
- CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
===========================================================================**/
+
#ifndef CLITKIMAGETOIMAGEGENERICFILTER_H
#define CLITKIMAGETOIMAGEGENERICFILTER_H
+
#include "clitkImageToImageGenericFilterBase.h"
namespace clitk {
// clitk
#include "clitkFilterBase.h"
+#include "clitkCropLikeImageFilter.h"
// itk
#include <itkPasteImageFilter.h>
void AddOrientationType(OrientationTypeEnumeration orientation);
void AddOrientationTypeString(std::string s);
void ClearOrientationType();
- void AddAngles(double a, double b);
+ void AddAnglesInRad(double a, double b);
+ void AddAnglesInDeg(double a, double b);
+ double GetAngle1InRad(int i) { return m_Angle1[i]; }
+ double GetAngle2InRad(int i) { return m_Angle2[i]; }
int GetNumberOfAngles();
std::string GetOrientationTypeString(int i) { return m_OrientationTypeString[i]; }
std::vector<std::string> & GetOrientationTypeString() { return m_OrientationTypeString; }
itkSetMacro(CombineWithOrFlag, bool);
itkBooleanMacro(CombineWithOrFlag);
+ itkGetConstMacro(FuzzyMapOnlyFlag, bool);
+ itkSetMacro(FuzzyMapOnlyFlag, bool);
+ itkBooleanMacro(FuzzyMapOnlyFlag);
+
+ typename FloatImageType::Pointer GetFuzzyMap() { return m_FuzzyMap; }
+
// I dont want to verify inputs information
virtual void VerifyInputInformation() { }
bool m_InverseOrientationFlag;
bool m_RemoveObjectFlag;
bool m_CombineWithOrFlag;
+ bool m_FuzzyMapOnlyFlag;
virtual void GenerateOutputInformation();
virtual void GenerateInputRequestedRegion();
typename ImageType::Pointer working_image;
typename ImageType::Pointer object_resampled;
typename FloatImageType::Pointer relPos;
+ typename FloatImageType::Pointer m_FuzzyMap;
ImagePointer input;
ImagePointer object;
#include "clitkAutoCropFilter.h"
#include "clitkResampleImageWithOptionsFilter.h"
#include "clitkBooleanOperatorLabelImageFilter.h"
+#include "clitkCropLikeImageFilter.h"
// itk
#include <deque>
CombineWithOrFlagOff();
VerboseStepFlagOff();
WriteStepFlagOff();
+ FuzzyMapOnlyFlagOff();
}
//--------------------------------------------------------------------
if (t == "NotSupTo") { AddOrientationType(SupTo); InverseOrientationFlagOn(); return; }
if (t == "NotInfTo") { AddOrientationType(InfTo); InverseOrientationFlagOn(); return; }
+ if (t == "Angle") return;
+
clitkExceptionMacro("Error, you must provide LeftTo,RightTo or AntTo,PostTo or SupTo,InfTo (or NotLeftTo, NotRightTo etc) but you give " << t);
}
//--------------------------------------------------------------------
template <class ImageType>
void
clitk::AddRelativePositionConstraintToLabelImageFilter<ImageType>::
-AddAngles(double a, double b)
+AddAnglesInRad(double a, double b)
{
- AddOrientationTypeString("Angle");
+ m_OrientationTypeString.push_back("Angle");
+ m_OrientationType.push_back(Angle);
m_Angle1.push_back(a);
m_Angle2.push_back(b);
}
//--------------------------------------------------------------------
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::AddRelativePositionConstraintToLabelImageFilter<ImageType>::
+AddAnglesInDeg(double a, double b)
+{
+ AddAnglesInRad(clitk::deg2rad(a), clitk::deg2rad(b));
+}
+//--------------------------------------------------------------------
+
+
//--------------------------------------------------------------------
template <class ImageType>
void
typedef itk::RelativePositionPropImageFilter<ImageType, FloatImageType> RelPosFilterType;
typename RelPosFilterType::Pointer relPosFilter;
- typename FloatImageType::Pointer m_FuzzyMap;
for(int i=0; i<GetNumberOfAngles(); i++) {
// Compute fuzzy map
relPosFilter = RelPosFilterType::New();
relPos = m_FuzzyMap;
StopCurrentStep<FloatImageType>(relPos);
+ if (GetFuzzyMapOnlyFlag()) return;
//--------------------------------------------------------------------
//--------------------------------------------------------------------
#include "clitkCommon.h"
#include "clitkBooleanOperatorLabelImageFilter.h"
-#include "clitkSegmentationUtils.h"
+#include "clitkBoundingBoxUtils.h"
namespace clitk {
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ======================================================================-====*/
+
+#ifndef CLITKBOUNDINGBOXUTILS_H
+#define CLITKBOUNDINGBOXUTILS_H
+
+// clitk
+#include "clitkCommon.h"
+
+// itk
+#include <itkBoundingBox.h>
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ template<class ImageType>
+ void ComputeBBFromImageRegion(const ImageType * image,
+ typename ImageType::RegionType region,
+ typename itk::BoundingBox<unsigned long,
+ ImageType::ImageDimension>::Pointer bb);
+
+ //--------------------------------------------------------------------
+ template<int Dimension>
+ void ComputeBBIntersection(typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbo,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi1,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi2);
+
+ //--------------------------------------------------------------------
+ template<int Dimension>
+ void ComputeBBUnion(typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbo,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi1,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi2);
+
+ //--------------------------------------------------------------------
+ template<class ImageType>
+ void ComputeRegionFromBB(const ImageType * image,
+ const typename itk::BoundingBox<unsigned long,
+ ImageType::ImageDimension>::Pointer bb,
+ typename ImageType::RegionType & region);
+
+} // end clitk namespace
+
+#include "clitkBoundingBoxUtils.txx"
+
+#endif
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ======================================================================-====*/
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ template<class ImageType>
+ void ComputeBBFromImageRegion(const ImageType * image,
+ typename ImageType::RegionType region,
+ typename itk::BoundingBox<unsigned long,
+ ImageType::ImageDimension>::Pointer bb) {
+ typedef typename ImageType::IndexType IndexType;
+ IndexType firstIndex;
+ IndexType lastIndex;
+ for(unsigned int i=0; i<image->GetImageDimension(); i++) {
+ firstIndex[i] = region.GetIndex()[i];
+ lastIndex[i] = firstIndex[i]+region.GetSize()[i];
+ }
+
+ typedef itk::BoundingBox<unsigned long,
+ ImageType::ImageDimension> BBType;
+ typedef typename BBType::PointType PointType;
+ PointType lastPoint;
+ PointType firstPoint;
+ image->TransformIndexToPhysicalPoint(firstIndex, firstPoint);
+ image->TransformIndexToPhysicalPoint(lastIndex, lastPoint);
+
+ bb->SetMaximum(lastPoint);
+ bb->SetMinimum(firstPoint);
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ template<int Dimension>
+ void ComputeBBIntersection(typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbo,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi1,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi2) {
+ typedef itk::BoundingBox<unsigned long, Dimension> BBType;
+ typedef typename BBType::PointType PointType;
+ PointType lastPoint;
+ PointType firstPoint;
+ for(unsigned int i=0; i<Dimension; i++) {
+ firstPoint[i] = std::max(bbi1->GetMinimum()[i], bbi2->GetMinimum()[i]);
+ lastPoint[i] = std::min(bbi1->GetMaximum()[i], bbi2->GetMaximum()[i]);
+ }
+ bbo->SetMaximum(lastPoint);
+ bbo->SetMinimum(firstPoint);
+ }
+ //--------------------------------------------------------------------
+
+
+ ///--------------------------------------------------------------------
+ template<int Dimension>
+ void ComputeBBUnion(typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbo,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi1,
+ typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi2) {
+ typedef itk::BoundingBox<unsigned long, Dimension> BBType;
+ typedef typename BBType::PointType PointType;
+ PointType lastPoint;
+ PointType firstPoint;
+ for(unsigned int i=0; i<Dimension; i++) {
+ firstPoint[i] = std::min(bbi1->GetMinimum()[i], bbi2->GetMinimum()[i]);
+ lastPoint[i] = std::max(bbi1->GetMaximum()[i], bbi2->GetMaximum()[i]);
+ }
+ bbo->SetMaximum(lastPoint);
+ bbo->SetMinimum(firstPoint);
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ template<class ImageType>
+ void ComputeRegionFromBB(const ImageType * image,
+ const typename itk::BoundingBox<unsigned long,
+ ImageType::ImageDimension>::Pointer bb,
+ typename ImageType::RegionType & region) {
+ // Types
+ typedef typename ImageType::IndexType IndexType;
+ typedef typename ImageType::PointType PointType;
+ typedef typename ImageType::RegionType RegionType;
+ typedef typename ImageType::SizeType SizeType;
+
+ // Region starting point
+ IndexType regionStart;
+ PointType start = bb->GetMinimum();
+ image->TransformPhysicalPointToIndex(start, regionStart);
+
+ // Region size
+ SizeType regionSize;
+ PointType maxs = bb->GetMaximum();
+ PointType mins = bb->GetMinimum();
+ for(unsigned int i=0; i<ImageType::ImageDimension; i++) {
+ regionSize[i] = lrint((maxs[i] - mins[i])/image->GetSpacing()[i]);
+ }
+
+ // Create region
+ region.SetIndex(regionStart);
+ region.SetSize(regionSize);
+ }
+ //--------------------------------------------------------------------
+
+
+} // end of namespace
+
#ifndef CLITKCROPLIKEIMAGEFILTER_H
#define CLITKCROPLIKEIMAGEFILTER_H
+// clitk
+#include "clitkBoundingBoxUtils.h"
+
+// itk
#include <itkImageToImageFilter.h>
namespace clitk {
}; // end class
//--------------------------------------------------------------------
+
+ //--------------------------------------------------------------------
+ // Convenient function
+ template<class ImageType>
+ typename ImageType::Pointer
+ ResizeImageLike(const ImageType * input,
+ const itk::ImageBase<ImageType::ImageDimension> * like,
+ typename ImageType::PixelType BG);
+
+ template<class ImageType>
+ typename ImageType::Pointer
+ ResizeImageLike(const ImageType * input,
+ typename itk::ImageBase<ImageType::ImageDimension>::RegionType * like,
+ typename ImageType::PixelType BG);
+
+ template<class ImageType>
+ typename ImageType::Pointer
+ ResizeImageLike(const ImageType * input,
+ typename itk::BoundingBox<unsigned long, ImageType::ImageDimension>::Pointer bb,
+ typename ImageType::PixelType BG);
+
+
} // end namespace clitk
//--------------------------------------------------------------------
#include "clitkCommon.h"
#include "clitkPasteImageFilter.h"
-// itk
-//#include "itkPasteImageFilter.h"
-
//--------------------------------------------------------------------
template <class ImageType>
clitk::CropLikeImageFilter<ImageType>::
}
//--------------------------------------------------------------------
-
-#endif //#define CLITKAUTOCROPFILTER
+
+//--------------------------------------------------------------------
+template<class ImageType>
+typename ImageType::Pointer
+clitk::ResizeImageLike(const ImageType * input,
+ const itk::ImageBase<ImageType::ImageDimension> * like,
+ typename ImageType::PixelType backgroundValue)
+{
+ typedef clitk::CropLikeImageFilter<ImageType> CropFilterType;
+ typename CropFilterType::Pointer cropFilter = CropFilterType::New();
+ cropFilter->SetInput(input);
+ cropFilter->SetCropLikeImage(like);
+ cropFilter->SetBackgroundValue(backgroundValue);
+ cropFilter->Update();
+ return cropFilter->GetOutput();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ImageType>
+typename ImageType::Pointer
+clitk::ResizeImageLike(const ImageType * input,
+ typename itk::ImageBase<ImageType::ImageDimension>::RegionType * region,
+ typename ImageType::PixelType backgroundValue)
+{
+ typename ImageType::Pointer output = ImageType::New();
+ output->CopyInformation(input);
+ output->SetRegions(region);
+ output->Allocate();
+ return clitk::ResizeImageLike<ImageType>(input, output, backgroundValue);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ImageType>
+typename ImageType::Pointer
+clitk::ResizeImageLike(const ImageType * input,
+ typename itk::BoundingBox<unsigned long, ImageType::ImageDimension>::Pointer bb,
+ typename ImageType::PixelType BG)
+{
+ typename ImageType::RegionType region;
+ clitk::ComputeRegionFromBB<ImageType>(input, bb, region);
+ typename ImageType::Pointer output = ImageType::New();
+ output->CopyInformation(input);
+ output->SetRegions(region);
+ output->Allocate();
+ return clitk::ResizeImageLike<ImageType>(input, output, BG);
+}
+//--------------------------------------------------------------------
+
+#endif //#define CLITKCROPLIKEIMAGEFILTER_TXX
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+#ifndef CLITKLABELIMAGEOVERLAPMEASUREFILTER_H
+#define CLITKRELATIVEPOSITIONANALYZERFILTER_H
+
+// clitk
+#include "clitkFilterBase.h"
+#include "clitkCropLikeImageFilter.h"
+#include "clitkSegmentationUtils.h"
+// #include "itkLabelOverlapMeasuresImageFilter.h"
+
+// itk
+#include <itkImageToImageFilter.h>
+#include <itkLabelStatisticsImageFilter.h>
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ /*
+ Analyze the relative position of a Target (mask image) according
+ to some Object, in a given Support. Indicate the main important
+ position of this Target according the Object.
+ */
+ //--------------------------------------------------------------------
+
+ template <class ImageType>
+ class LabelImageOverlapMeasureFilter:
+ public virtual FilterBase,
+ public itk::ImageToImageFilter<ImageType, ImageType>
+ {
+
+ public:
+ /** Standard class typedefs. */
+ typedef itk::ImageToImageFilter<ImageType, ImageType> Superclass;
+ typedef LabelImageOverlapMeasureFilter<ImageType> Self;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ /** Method for creation through the object factory. */
+ itkNewMacro(Self);
+
+ /** Run-time type information (and related methods). */
+ itkTypeMacro(LabelImageOverlapMeasureFilter, ImageToImageFilter);
+
+ /** Some convenient typedefs. */
+ typedef typename ImageType::ConstPointer ImageConstPointer;
+ typedef typename ImageType::Pointer ImagePointer;
+ typedef typename ImageType::RegionType RegionType;
+ typedef typename ImageType::PixelType PixelType;
+ typedef typename ImageType::SpacingType SpacingType;
+ typedef typename ImageType::SizeType SizeType;
+ typedef typename ImageType::IndexType IndexType;
+ typedef typename ImageType::PointType PointType;
+
+ /** ImageDimension constants */
+ itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension);
+ FILTERBASE_INIT;
+
+ // Options
+ itkGetConstMacro(BackgroundValue, PixelType);
+ itkSetMacro(BackgroundValue, PixelType);
+
+ itkGetConstMacro(Label1, PixelType);
+ itkSetMacro(Label1, PixelType);
+
+ itkGetConstMacro(Label2, PixelType);
+ itkSetMacro(Label2, PixelType);
+
+ // I dont want to verify inputs information
+ virtual void VerifyInputInformation() { }
+
+ protected:
+ LabelImageOverlapMeasureFilter();
+ virtual ~LabelImageOverlapMeasureFilter() {}
+
+ PixelType m_BackgroundValue;
+ PixelType m_Label1;
+ PixelType m_Label2;
+ ImagePointer m_Input1;
+ ImagePointer m_Input2;
+
+ virtual void GenerateOutputInformation();
+ virtual void GenerateInputRequestedRegion();
+ virtual void GenerateData();
+
+ private:
+ LabelImageOverlapMeasureFilter(const Self&); //purposely not implemented
+ void operator=(const Self&); //purposely not implemented
+
+ }; // end class
+ //--------------------------------------------------------------------
+
+} // end namespace clitk
+//--------------------------------------------------------------------
+
+#include "clitkLabelImageOverlapMeasureFilter.txx"
+
+#endif
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+//--------------------------------------------------------------------
+template <class ImageType>
+clitk::LabelImageOverlapMeasureFilter<ImageType>::
+LabelImageOverlapMeasureFilter():
+ clitk::FilterBase(),
+ itk::ImageToImageFilter<ImageType, ImageType>()
+{
+ this->SetNumberOfRequiredInputs(2);
+ SetLabel1(1);
+ SetLabel2(1);
+ SetBackgroundValue(0);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::LabelImageOverlapMeasureFilter<ImageType>::
+GenerateOutputInformation()
+{
+ // DD("GenerateOutputInformation");
+ //ImagePointer input = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ // ImagePointer outputImage = this->GetOutput(0);
+ // outputImage->SetRegions(outputImage->GetLargestPossibleRegion());
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::LabelImageOverlapMeasureFilter<ImageType>::
+GenerateInputRequestedRegion()
+{
+ // DD("GenerateInputRequestedRegion");
+ // Call default
+ itk::ImageToImageFilter<ImageType, ImageType>::GenerateInputRequestedRegion();
+ // Get input pointers and set requested region to common region
+ ImagePointer input1 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ ImagePointer input2 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(1));
+ input1->SetRequestedRegion(input1->GetLargestPossibleRegion());
+ input2->SetRequestedRegion(input2->GetLargestPossibleRegion());
+}
+//--------------------------------------------------------------------
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::LabelImageOverlapMeasureFilter<ImageType>::
+GenerateData()
+{
+ // DD("GenerateData");
+
+ // Get input pointer
+ m_Input1 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ m_Input2 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(1));
+ static const unsigned int dim = ImageType::ImageDimension;
+
+ // Compute union of bounding boxes
+ typedef itk::BoundingBox<unsigned long, dim> BBType;
+ typename BBType::Pointer bb1 = BBType::New();
+ ComputeBBFromImageRegion<ImageType>(m_Input1, m_Input1->GetLargestPossibleRegion(), bb1);
+ typename BBType::Pointer bb2 = BBType::New();
+ ComputeBBFromImageRegion<ImageType>(m_Input2, m_Input2->GetLargestPossibleRegion(), bb2);
+ typename BBType::Pointer bbo = BBType::New();
+ ComputeBBUnion<dim>(bbo, bb1, bb2);
+
+ // Resize like the union
+ ImagePointer input1 = clitk::ResizeImageLike<ImageType>(m_Input1, bbo, GetBackgroundValue());
+ ImagePointer input2 = clitk::ResizeImageLike<ImageType>(m_Input2, bbo, GetBackgroundValue());
+
+ // Compute overlap image
+ ImagePointer image_union = clitk::Clone<ImageType>(input1);
+ ImagePointer image_intersection = clitk::Clone<ImageType>(input1);
+ clitk::Or<ImageType>(image_union, input2, GetBackgroundValue());
+ clitk::And<ImageType>(image_intersection, input2, GetBackgroundValue());
+
+ writeImage<ImageType>(image_union, "union.mha");
+ writeImage<ImageType>(image_intersection, "intersection.mha");
+
+ // Compute size
+ typedef itk::LabelStatisticsImageFilter<ImageType, ImageType> StatFilterType;
+ typename StatFilterType::Pointer statFilter = StatFilterType::New();
+ statFilter->SetInput(image_union);
+ statFilter->SetLabelInput(image_union);
+ statFilter->Update();
+ int u = statFilter->GetCount(GetLabel1());
+
+ statFilter->SetInput(image_intersection);
+ statFilter->SetLabelInput(image_intersection);
+ statFilter->Update();
+ int inter = statFilter->GetCount(GetLabel1());
+
+ statFilter->SetInput(m_Input1);
+ statFilter->SetLabelInput(m_Input1);
+ statFilter->Update();
+ int in1 = statFilter->GetCount(GetLabel1());
+
+ statFilter->SetInput(m_Input2);
+ statFilter->SetLabelInput(m_Input2);
+ statFilter->Update();
+ int in2 = statFilter->GetCount(GetLabel1());
+
+ std::cout << in1 << " " << in2 << " " << inter << " " << u << " " << (double)inter/(double)u << std::endl;
+}
+//--------------------------------------------------------------------
+
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+#ifndef CLITKRELATIVEPOSITIONANALYZERFILTER_H
+#define CLITKRELATIVEPOSITIONANALYZERFILTER_H
+
+// clitk
+#include "clitkFilterWithAnatomicalFeatureDatabaseManagement.h"
+#include "clitkFilterBase.h"
+#include "clitkSliceBySliceRelativePositionFilter.h"
+#include "clitkRelativePositionDataBase.h"
+
+// itk
+#include <itkImageToImageFilter.h>
+#include <itkLabelStatisticsImageFilter.h>
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ /*
+ Analyze the relative position of a Target (mask image) according
+ to some Object (mask image), in a given Support (mask
+ image). Compute the optimal threshold allowing to remove the
+ maximal area from the Support without removing area belonging to
+ the Target.
+ */
+ //--------------------------------------------------------------------
+
+ template <class ImageType>
+ class RelativePositionAnalyzerFilter:
+ public itk::ImageToImageFilter<ImageType, ImageType>
+ {
+
+ public:
+ /** Standard class typedefs. */
+ typedef itk::ImageToImageFilter<ImageType, ImageType> Superclass;
+ typedef RelativePositionAnalyzerFilter<ImageType> Self;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ /** Method for creation through the object factory. */
+ itkNewMacro(Self);
+
+ /** Run-time type information (and related methods). */
+ itkTypeMacro(RelativePositionAnalyzerFilter, ImageToImageFilter);
+
+ /** Some convenient typedefs. */
+ typedef typename ImageType::ConstPointer ImageConstPointer;
+ typedef typename ImageType::Pointer ImagePointer;
+ typedef typename ImageType::RegionType RegionType;
+ typedef typename ImageType::PixelType PixelType;
+ typedef typename ImageType::SpacingType SpacingType;
+ typedef typename ImageType::SizeType SizeType;
+ typedef typename ImageType::IndexType IndexType;
+ typedef typename ImageType::PointType PointType;
+
+ /** ImageDimension constants */
+ itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension);
+ FILTERBASE_INIT;
+ typedef itk::Image<float, ImageDimension> FloatImageType;
+
+ /** Input : initial image and object */
+ void SetInputSupport(const ImageType * image);
+ void SetInputObject(const ImageType * image);
+ void SetInputTarget(const ImageType * image);
+
+ // Input
+ // supportname, objectname multiple targetname
+
+ // Options
+ itkGetConstMacro(BackgroundValue, PixelType);
+ itkSetMacro(BackgroundValue, PixelType);
+
+ itkGetConstMacro(ForegroundValue, PixelType);
+ itkSetMacro(ForegroundValue, PixelType);
+
+ clitk::RelativePositionDirectionType & GetDirection() { return m_Direction; }
+ void SetDirection(clitk::RelativePositionDirectionType & d) { m_Direction = d; }
+
+ itkGetConstMacro(NumberOfBins, int);
+ itkSetMacro(NumberOfBins, int);
+
+ itkGetConstMacro(AreaLossTolerance, double);
+ itkSetMacro(AreaLossTolerance, double);
+
+ itkGetConstMacro(SupportSize, int);
+ itkGetConstMacro(TargetSize, int);
+ itkGetConstMacro(SizeWithThreshold, int);
+ itkGetConstMacro(SizeWithReverseThreshold, int);
+
+ itkGetConstMacro(Info, clitk::RelativePositionInformationType);
+ itkGetConstMacro(InfoReverse, clitk::RelativePositionInformationType);
+
+ // For debug
+ void PrintOptions();
+
+ // Print output
+ void Print(std::ostream & os=std::cout);
+
+ // I dont want to verify inputs information
+ virtual void VerifyInputInformation() { }
+
+ protected:
+ RelativePositionAnalyzerFilter();
+ virtual ~RelativePositionAnalyzerFilter() {}
+
+ itkSetMacro(SupportSize, int);
+ itkSetMacro(TargetSize, int);
+ itkSetMacro(SizeWithThreshold, int);
+ itkSetMacro(SizeWithReverseThreshold, int);
+
+ PixelType m_BackgroundValue;
+ PixelType m_ForegroundValue;
+ ImagePointer m_Support;
+ ImagePointer m_Object;
+ ImagePointer m_Target;
+ int m_NumberOfBins;
+ double m_AreaLossTolerance;
+ int m_SupportSize;
+ int m_TargetSize;
+ int m_SizeWithReverseThreshold;
+ int m_SizeWithThreshold;
+ clitk::RelativePositionDirectionType m_Direction;
+ clitk::RelativePositionInformationType m_Info;
+ clitk::RelativePositionInformationType m_InfoReverse;
+
+ virtual void GenerateOutputInformation();
+ virtual void GenerateData();
+
+ typename FloatImageType::Pointer
+ ComputeFuzzyMap(ImageType * object, ImageType * target, ImageType * support, double angle);
+
+ void
+ ComputeOptimalThresholds(FloatImageType * map, ImageType * target, int bins, double tolerance,
+ double & threshold, double & reverseThreshold);
+ private:
+ RelativePositionAnalyzerFilter(const Self&); //purposely not implemented
+ void operator=(const Self&); //purposely not implemented
+
+ }; // end class
+ //--------------------------------------------------------------------
+
+} // end namespace clitk
+//--------------------------------------------------------------------
+
+#include "clitkRelativePositionAnalyzerFilter.txx"
+
+#endif
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+//--------------------------------------------------------------------
+template <class ImageType>
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+RelativePositionAnalyzerFilter():
+ itk::ImageToImageFilter<ImageType, ImageType>()
+{
+ this->SetNumberOfRequiredInputs(3); // Input : support, object, target
+ SetBackgroundValue(0);
+ SetForegroundValue(1);
+ SetNumberOfBins(100);
+ SetAreaLossTolerance(0.01);
+ SetSupportSize(0);
+ SetTargetSize(0);
+ SetSizeWithThreshold(0);
+ SetSizeWithReverseThreshold(0);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+SetInputSupport(const ImageType * image)
+{
+ // Process object is not const-correct so the const casting is required.
+ this->SetNthInput(0, const_cast<ImageType *>(image));
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+SetInputObject(const ImageType * image)
+{
+ // Process object is not const-correct so the const casting is required.
+ this->SetNthInput(1, const_cast<ImageType *>(image));
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+SetInputTarget(const ImageType * image)
+{
+ // Process object is not const-correct so the const casting is required.
+ this->SetNthInput(2, const_cast<ImageType *>(image));
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+PrintOptions()
+{
+ DD("TODO");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+GenerateOutputInformation()
+{
+ ImagePointer input = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ ImagePointer outputImage = this->GetOutput(0);
+ outputImage->SetRegions(outputImage->GetLargestPossibleRegion());
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+GenerateData()
+{
+ static const unsigned int dim = ImageType::ImageDimension;
+
+ ImagePointer temp = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ m_Object = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(1));
+ m_Target = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(2));
+
+ // Remove object from support (keep initial image)
+ m_Support = clitk::Clone<ImageType>(temp);
+ clitk::AndNot<ImageType>(m_Support, m_Object, GetBackgroundValue());
+
+ // Define filter to compute statics on mask image
+ typedef itk::LabelStatisticsImageFilter<ImageType, ImageType> StatFilterType;
+ typename StatFilterType::Pointer statFilter = StatFilterType::New();
+
+ // Compute the initial support size
+ statFilter->SetInput(m_Support);
+ statFilter->SetLabelInput(m_Support);
+ statFilter->Update();
+ SetSupportSize(statFilter->GetCount(GetForegroundValue()));
+ // DD(GetSupportSize());
+
+ // Compute the initial target size
+ ImagePointer s = clitk::ResizeImageLike<ImageType>(m_Support, m_Target, GetBackgroundValue());
+ statFilter->SetInput(s);
+ statFilter->SetLabelInput(m_Target);
+ statFilter->Update();
+ SetTargetSize(statFilter->GetCount(GetForegroundValue()));
+ // DD(GetTargetSize());
+
+ //
+ int bins = GetNumberOfBins();
+ double tolerance = GetAreaLossTolerance();
+
+ // Compute Fuzzy map
+ double angle = GetDirection().angle1;
+ typename FloatImageType::Pointer map = ComputeFuzzyMap(m_Object, m_Target, m_Support, angle);
+ writeImage<FloatImageType>(map, "fuzzy_"+toString(clitk::rad2deg(angle))+".mha");
+
+ // Compute the optimal thresholds (direct and inverse)
+ double mThreshold=0.0;
+ double mReverseThreshold=1.0;
+ ComputeOptimalThresholds(map, m_Target, bins, tolerance, mThreshold, mReverseThreshold);
+
+ // Use the threshold to compute new support
+ int s1 = GetSupportSize();
+ if (mThreshold > 0.0) {
+ ImagePointer support1 =
+ clitk::SliceBySliceRelativePosition<ImageType>(m_Support, m_Object, 2,
+ mThreshold,
+ angle,false, // inverseFlag
+ false, // uniqueConnectedComponent
+ -1, true,
+ false);//singleObjectCCL
+ // Compute the new support size
+ statFilter->SetInput(support1);
+ statFilter->SetLabelInput(support1);
+ statFilter->Update();
+ s1 = statFilter->GetCount(GetForegroundValue());
+ }
+
+ int s2 = GetSupportSize();
+ if (mReverseThreshold < 1.0) {
+ ImagePointer support2 =
+ clitk::SliceBySliceRelativePosition<ImageType>(m_Support, m_Object, 2,
+ mReverseThreshold,
+ angle,true,// inverseFlag
+ false, // uniqueConnectedComponent
+ -1, true,
+ false); //singleObjectCCL
+ // Compute the new support size
+ statFilter = StatFilterType::New();
+ statFilter->SetInput(support2);
+ statFilter->SetLabelInput(support2);
+ statFilter->Update();
+ s2 = statFilter->GetCount(GetForegroundValue());
+ }
+
+ // Set results values
+ m_Info.threshold = mThreshold;
+ m_Info.sizeAfterThreshold = s1;
+ m_Info.sizeBeforeThreshold = GetSupportSize();
+ m_Info.sizeReference = GetTargetSize();
+ m_InfoReverse.threshold = mReverseThreshold;
+ m_InfoReverse.sizeAfterThreshold = s2;
+ m_InfoReverse.sizeBeforeThreshold = GetSupportSize();
+ m_InfoReverse.sizeReference = GetTargetSize();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+typename clitk::RelativePositionAnalyzerFilter<ImageType>::FloatImageType::Pointer
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+ComputeFuzzyMap(ImageType * object, ImageType * target, ImageType * support, double angle)
+{
+ typedef clitk::SliceBySliceRelativePositionFilter<ImageType> SliceRelPosFilterType;
+ typedef typename SliceRelPosFilterType::FloatImageType FloatImageType;
+ typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New();
+ sliceRelPosFilter->VerboseStepFlagOff();
+ sliceRelPosFilter->WriteStepFlagOff();
+ sliceRelPosFilter->SetInput(support);
+ sliceRelPosFilter->SetInputObject(object);
+ sliceRelPosFilter->SetDirection(2);
+ sliceRelPosFilter->SetIntermediateSpacingFlag(false);
+ //sliceRelPosFilter->AddOrientationTypeString(orientation);
+ sliceRelPosFilter->AddAnglesInRad(angle, 0.0);
+ sliceRelPosFilter->FuzzyMapOnlyFlagOn(); // do not threshold, only compute the fuzzy map
+ // sliceRelPosFilter->PrintOptions();
+ sliceRelPosFilter->Update();
+ typename FloatImageType::Pointer map = sliceRelPosFilter->GetFuzzyMap();
+
+ // Resize map like object to allow SetBackground
+ map = clitk::ResizeImageLike<FloatImageType>(map, object, GetBackgroundValue());
+
+ // Remove initial object from the fuzzy map
+ map = clitk::SetBackground<FloatImageType, ImageType>(map, object, GetForegroundValue(), 0.0, true);
+
+ // Resize the fuzzy map like the target, put 2.0 when outside
+ map = clitk::ResizeImageLike<FloatImageType>(map, target, 2.0); // Put 2.0 when out of initial map
+
+ // end
+ return map;
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionAnalyzerFilter<ImageType>::
+ComputeOptimalThresholds(FloatImageType * map, ImageType * target, int bins, double tolerance,
+ double & threshold, double & reverseThreshold)
+{
+ // Get the histogram of fuzzy values inside the target image
+ typedef itk::LabelStatisticsImageFilter<FloatImageType, ImageType> FloatStatFilterType;
+ typename FloatStatFilterType::Pointer f = FloatStatFilterType::New();
+ f->SetInput(map);
+ f->SetLabelInput(target);
+ f->UseHistogramsOn();
+ f->SetHistogramParameters(bins, 0.0, 1.1);
+ f->Update();
+ int count = f->GetCount(GetForegroundValue());
+ typename FloatStatFilterType::HistogramPointer h = f->GetHistogram(GetForegroundValue());
+
+ // Debug : dump histogram
+ static int i=0;
+ std::ofstream histogramFile(std::string("fuzzy_histo_"+toString(i)+".txt").c_str());
+ for(int j=0; j<bins; j++) {
+ histogramFile << h->GetMeasurement(j,0)
+ << "\t" << h->GetFrequency(j)
+ << "\t" << (double)h->GetFrequency(j)/(double)count << std::endl;
+ }
+ histogramFile.close();
+ i++;
+
+ // Analyze the histogram (direct)
+ double sum = 0.0;
+ bool found = false;
+ threshold = 0.0;
+ for(int j=0; j<bins; j++) {
+ sum += ((double)h->GetFrequency(j)/(double)count);
+ if ((!found) && (sum > tolerance)) {
+ if (j==0) threshold = h->GetBinMin(0,j);
+ else threshold = h->GetBinMin(0,j-1); // the last before reaching the threshold
+ found = true;
+ }
+ }
+
+ // Analyze the histogram (reverse)
+ sum = 0.0;
+ found = false;
+ reverseThreshold = 1.0;
+ for(int j=bins-1; j>=0; j--) {
+ sum += ((double)h->GetFrequency(j)/(double)count);
+ if ((!found) && (sum > tolerance)) {
+ if (j==bins-1) reverseThreshold = h->GetBinMax(0,j);
+ else reverseThreshold = h->GetBinMax(0,j+1);
+ found = true;
+ }
+ }
+}
+//--------------------------------------------------------------------
+
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+#ifndef CLITKRELATIVEPOSITIONDATABASE_CXX
+#define CLITKRELATIVEPOSITIONDATABASE_CXX
+
+// clitk
+#include "clitkRelativePositionDataBase.h"
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ void RelativePositionDataBase::ReadIndex(std::istream & is, IndexType & index)
+ {
+ is >> index.patient;
+ is >> index.station;
+ is >> index.object;
+ is >> index.direction.angle1;
+ index.direction.angle1 = clitk::deg2rad(index.direction.angle1);
+ is >> index.direction.angle2;
+ index.direction.angle2 = clitk::deg2rad(index.direction.angle2);
+ std::string s;
+ is >> s;
+ if (s=="true") index.direction.notFlag = true;
+ else index.direction.notFlag = false;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ void RelativePositionDataBase::ReadInformation(std::istream & is, RelativePositionInformationType & v)
+ {
+ is >> v.threshold;
+ is >> v.sizeBeforeThreshold;
+ is >> v.sizeAfterThreshold;
+ is >> v.sizeReference;
+ }
+ //--------------------------------------------------------------------
+
+ //--------------------------------------------------------------------
+ void RelativePositionDataBase::Read(const std::string & filename)
+ {
+ std::ifstream is;
+ openFileForReading(is, filename);
+
+ std::string s;
+ IndexType index;
+ RelativePositionInformationType v;
+ while (is) {
+ skipComment(is); /* FIXME Read Index etc */
+ ReadIndex(is, index);
+ ReadInformation(is, v);
+
+ if (is) {// FIXME INSERT
+ // Set in station
+ if (m_DB.find(index.station) == m_DB.end()) {
+ MapByObjectType s;
+ m_DB[index.station] = s;
+ }
+ MapByObjectType & s = m_DB[index.station];
+
+ // Get Direction map from Object
+ if (s.find(index.object) == s.end()) {
+ MapByDirectionType r;
+ s[index.object] = r;
+ }
+ MapByDirectionType & r = s[index.object];
+
+ // Get Patient map from Direction
+ if (r.find(index.direction) == r.end()) {
+ MapByPatientType q;
+ r[index.direction] = q;
+ }
+ MapByPatientType & q = r[index.direction];
+
+ // Set value by patient
+ q[index.patient] = v;
+
+ // Debug
+ // index.Println();
+ // GetMapByPatient(index).find(index.patient)->second.Println();
+
+ } // End insertion
+ } // end loop reading
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ const RelativePositionDataBase::MapByDirectionType &
+ RelativePositionDataBase::GetMapByDirection(const IndexType & index) const
+ {
+ const MapByObjectType & a = GetMapByObject(index.station);
+ if (a.find(index.object) == a.end()) {
+ clitkExceptionMacro("Could not find index in DB (object= '" << index.object << "' not found).");
+ }
+ return a.find(index.object)->second;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ const RelativePositionDataBase::MapByObjectType &
+ RelativePositionDataBase::GetMapByObject(const std::string & station) const
+ {
+ if (m_DB.find(station) == m_DB.end()) {
+ clitkExceptionMacro("Could not find index in DB (station= '" << station << "' not found).");
+ }
+ return m_DB.find(station)->second;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ const RelativePositionDataBase::MapByPatientType &
+ RelativePositionDataBase::GetMapByPatient(const IndexType & index) const
+ {
+ const MapByDirectionType & a = GetMapByDirection(index);
+ if (a.find(index.direction) == a.end()) {
+ std::ostringstream s;
+ index.direction.Print(s);
+ clitkExceptionMacro("Could not find index in DB (direction= '" << s.str() << "' not found).");
+ }
+ return a.find(index.direction)->second;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ const RelativePositionInformationType &
+ RelativePositionDataBase::GetInformation(const IndexType & index) const
+ {
+ const RelativePositionDataBase::MapByPatientType & a = GetMapByPatient(index);
+ if (a.find(index.patient) == a.end()) {
+ clitkExceptionMacro("Could not find index in DB (patient= '" << index.patient << "' not found).");
+ }
+ return a.find(index.patient)->second;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ int RelativePositionDataBase::GetNumberOfPatient(const IndexType & index) const
+ {
+ const MapByPatientType & o = GetMapByPatient(index);
+ return o.size();
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ std::vector<std::string> & RelativePositionDataBase::GetListOfPatients(const IndexType & index) const
+ {
+ const MapByPatientType & o = GetMapByPatient(index);
+ MapByPatientType::const_iterator iter = o.begin();
+ std::vector<std::string> * v = new std::vector<std::string>;
+ MapToVecFirst(o, *v);
+ return *v;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ double RelativePositionDataBase::GetAreaGain(const IndexType & index) const
+ {
+ // FIXME change name
+ const RelativePositionInformationType & v = GetInformation(index);
+ return v.sizeAfterThreshold/v.sizeBeforeThreshold;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ double RelativePositionDataBase::GetThreshold(const IndexType & index) const
+ {
+ const RelativePositionInformationType & v = GetInformation(index);
+ return v.threshold;
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ void
+ RelativePositionDataBase::GetListOfObjects(const std::string & station, std::vector<std::string> & objects) const
+ {
+ const MapByObjectType & a = GetMapByObject(station);
+ MapToVecFirst(a, objects);
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ void
+ RelativePositionDataBase::GetListOfDirections(const std::string & station,
+ const std::string & object,
+ std::vector<RelativePositionDirectionType> & directions) const
+ {
+ IndexType i;
+ i.station = station;
+ i.object = object;
+ const MapByDirectionType & n = GetMapByDirection(i);
+ MapToVecFirst(n, directions);
+ }
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ bool RelativePositionDataBase::CheckIndex(const IndexType & index) const
+ {
+ try {
+ const RelativePositionInformationType & m = GetInformation(index);
+ } catch (clitk::ExceptionObject e) {
+ // std::cout << e.what() << std::endl;
+ return false;
+ }
+ return true;
+ }
+ //--------------------------------------------------------------------
+
+
+} // end namespace clitk
+//--------------------------------------------------------------------
+
+#endif
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+#ifndef CLITKRELATIVEPOSITIONDATABASE_H
+#define CLITKRELATIVEPOSITIONDATABASE_H
+
+// clitk
+#include "clitkCommon.h"
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ /*
+ FIXME
+ */
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ class RelativePositionDirectionType {
+ public:
+ double angle1;
+ double angle2;
+ bool notFlag;
+
+ void Print(std::ostream & os = std::cout) const {
+ os << clitk::rad2deg(angle1) << " " << clitk::rad2deg(angle2) << " ";
+ if (notFlag) os << "true";
+ else os << "false";
+ os << " ";
+ }
+
+ void PrintOptions(std::ostream & os = std::cout) const {
+ os << "angle1 = " << clitk::rad2deg(angle1) << std::endl
+ << "angle2 = " << clitk::rad2deg(angle2) << std::endl;
+ if (notFlag) os << "inverse" << std::endl;
+ }
+
+ void Println(std::ostream & os = std::cout) const {
+ Print(os);
+ os << std::endl;
+ }
+
+ bool operator< (const RelativePositionDirectionType &compare) const
+ {
+ if (angle1 < compare.angle1) return true;
+ if (angle1 > compare.angle1) return false;
+
+ if (angle2 < compare.angle2) return true;
+ if (angle2 > compare.angle2) return false;
+
+ if (notFlag == true) {
+ if (compare.notFlag == false) return true;
+ else return false;
+ }
+ return false;
+ }
+ };
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ class RelativePositionDataBaseIndexType {
+ public:
+ std::string patient;
+ std::string station;
+ std::string object;
+ RelativePositionDirectionType direction;
+ void Print(std::ostream & os = std::cout) const {
+ os << patient << " " << station << " " << object << " ";
+ direction.Print(os);
+ }
+ void Println(std::ostream & os = std::cout) const {
+ Print(os);
+ os << std::endl;
+ }
+ };
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ class RelativePositionInformationType {
+ public:
+ double threshold;
+ int sizeBeforeThreshold;
+ int sizeAfterThreshold;
+ int sizeReference;
+ void Print(std::ostream & os = std::cout) const {
+ os << threshold << " " << sizeBeforeThreshold << " "
+ << sizeAfterThreshold << " " << sizeReference;
+ }
+ void Println(std::ostream & os = std::cout) const {
+ Print(os);
+ os << std::endl;
+ }
+ };
+ //--------------------------------------------------------------------
+
+
+ //--------------------------------------------------------------------
+ class RelativePositionDataBase {
+
+ public:
+ RelativePositionDataBase() {}
+ ~RelativePositionDataBase() {}
+
+ typedef RelativePositionDataBaseIndexType IndexType;
+
+ void Read(const std::string & filename);
+ double GetAreaGain(const IndexType & index) const;
+ double GetThreshold(const IndexType & index) const;
+ int GetNumberOfPatient(const IndexType & index) const;
+ std::vector<std::string> & GetListOfPatients(const IndexType & index) const;
+ void GetListOfObjects(const std::string & station, std::vector<std::string> & objects) const;
+ void GetListOfDirections(const std::string & station,
+ const std::string & object,
+ std::vector<RelativePositionDirectionType> & directions) const;
+ bool CheckIndex(const IndexType & index) const;
+
+ protected:
+ typedef std::map<std::string, RelativePositionInformationType> MapByPatientType;
+ typedef std::map<RelativePositionDirectionType, MapByPatientType> MapByDirectionType;
+ typedef std::map<std::string, MapByDirectionType> MapByObjectType;
+ typedef std::map<std::string, MapByObjectType> MapByStationType;
+ MapByStationType m_DB;
+
+ void ReadIndex(std::istream & is, IndexType & index);
+ void ReadInformation(std::istream & is, RelativePositionInformationType & v);
+
+ const MapByDirectionType & GetMapByDirection(const IndexType & index) const;
+ const MapByPatientType & GetMapByPatient(const IndexType & index) const;
+ const RelativePositionInformationType & GetInformation(const IndexType & index) const;
+ const MapByObjectType & GetMapByObject(const std::string & station) const;
+
+ }; // end class
+ //--------------------------------------------------------------------
+
+} // end namespace clitk
+//--------------------------------------------------------------------
+
+#endif
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+#ifndef CLITKRELATIVEPOSITIONANALYZERFILTER_H
+#define CLITKRELATIVEPOSITIONANALYZERFILTER_H
+
+// clitk
+#include "clitkFilterWithAnatomicalFeatureDatabaseManagement.h"
+#include "clitkFilterBase.h"
+#include "clitkRelativePositionDataBase.h"
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ /*
+ Load a database of relative positions. Analyze it and provide
+ common relative position for each patient.
+ */
+ //--------------------------------------------------------------------
+
+ class RelativePositionDataBaseAnalyzerFilter:
+ public virtual clitk::FilterBase,
+ public clitk::FilterWithAnatomicalFeatureDatabaseManagement
+ {
+
+ public:
+ RelativePositionDataBaseAnalyzerFilter();
+ virtual ~RelativePositionDataBaseAnalyzerFilter() {}
+
+ // Input
+ itkGetConstMacro(DatabaseFilename, std::string);
+ itkSetMacro(DatabaseFilename, std::string);
+ itkGetConstMacro(StationName, std::string);
+ itkSetMacro(StationName, std::string);
+ itkGetConstMacro(ObjectName, std::string);
+ itkSetMacro(ObjectName, std::string);
+ itkGetConstMacro(OutputFilename, std::string);
+ itkSetMacro(OutputFilename, std::string);
+
+ // For debug
+ void PrintOptions();
+
+ // Go !
+ virtual void Update();
+
+ protected:
+ std::string m_DatabaseFilename;
+ std::string m_StationName;
+ std::string m_ObjectName;
+ std::string m_OutputFilename;
+ clitk::RelativePositionDataBase db;
+
+ bool ComputeOptimalThreshold(RelativePositionDataBaseIndexType & index, double & threshold);
+
+ }; // end class
+ //--------------------------------------------------------------------
+
+} // end namespace clitk
+//--------------------------------------------------------------------
+
+#include "clitkRelativePositionDataBaseAnalyzerFilter.txx"
+
+#endif
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+//--------------------------------------------------------------------
+clitk::RelativePositionDataBaseAnalyzerFilter::
+RelativePositionDataBaseAnalyzerFilter():
+ clitk::FilterBase(),
+ clitk::FilterWithAnatomicalFeatureDatabaseManagement()
+{
+ VerboseFlagOff();
+ SetDatabaseFilename("default.dbrp");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+void
+clitk::RelativePositionDataBaseAnalyzerFilter::
+PrintOptions()
+{
+ DD("TODO");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+void
+clitk::RelativePositionDataBaseAnalyzerFilter::
+Update()
+{
+ // Load DB of relative position
+ db.Read(GetDatabaseFilename());
+
+ // Get list of objects, list of orientation
+ std::vector<std::string> m_ListOfObjects;
+ db.GetListOfObjects(GetStationName(), m_ListOfObjects);
+
+ // Set initial index in the DB
+ clitk::RelativePositionDataBase::IndexType index;
+ index.station = GetStationName();
+
+ // Loop over objects
+ std::vector<double> m_ListOfThresholds;
+ for(int i=0; i<m_ListOfObjects.size(); i++) {
+ // DD(i);
+ // DD(m_ListOfObjects[i]);
+ // Set current index
+ index.object = m_ListOfObjects[i];
+ // Get the list of direction
+ std::vector<clitk::RelativePositionDirectionType> m_ListOfDirections;
+ db.GetListOfDirections(GetStationName(), index.object, m_ListOfDirections);
+
+ // Loop over direction
+ for(int j=0; j<m_ListOfDirections.size(); j++) {
+ // DD(j);
+ // m_ListOfDirections[j].Println();
+ // Set current index
+ index.direction = m_ListOfDirections[j];
+ // Compute the best RelPos parameters
+ double threshold;
+ bool ok = ComputeOptimalThreshold(index, threshold);
+ m_ListOfThresholds.push_back(threshold);
+
+ // Print debug FIXME
+ if (ok) {
+ /*std::cout << m_ListOfObjects[i] << " ";
+ m_ListOfDirections[j].Print();
+ std::cout << " " << threshold << " " << ok << std::endl;
+ */
+ std::cout << "# -----------------------" << std::endl
+ << "object = " << m_ListOfObjects[i] << std::endl;
+ m_ListOfDirections[j].PrintOptions();
+ std::cout << "threshold = " << threshold << std::endl
+ << "sliceBySlice" << std::endl << std::endl; // FIXME spacing ?
+ }
+ }
+ }
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+bool
+clitk::RelativePositionDataBaseAnalyzerFilter::
+ComputeOptimalThreshold(RelativePositionDataBaseIndexType & index, double & threshold)
+{
+ // Get list of patient
+ std::vector<std::string> & ListOfPatients = db.GetListOfPatients(index);
+ // DD(ListOfPatients.size());
+ // index.direction.Println();
+
+ // For a given station, object, direction
+ bool stop=false;
+ int i=0;
+ if (index.direction.notFlag) threshold = 0.0;
+ else threshold = 1.0;
+ while (!stop && (i<ListOfPatients.size())) {
+ index.patient = ListOfPatients[i];
+ // Check index
+ if (!db.CheckIndex(index)) {
+ std::cout << "Warning index does not exist in the DB. index = ";
+ index.Println(std::cout);
+ }
+ else {
+ if (index.direction.notFlag) threshold = std::max(db.GetThreshold(index), threshold);
+ else threshold = std::min(db.GetThreshold(index), threshold);
+ }
+ ++i;
+ } // end while
+
+ if (index.direction.notFlag) {
+ if (threshold >=1) return false; // not useful
+ }
+ else {
+ if (threshold <=0) return false; // not useful
+ }
+ return true;
+}
+//--------------------------------------------------------------------
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+#ifndef CLITKRelativePositionDataBaseBuilderFILTER_H
+#define CLITKRelativePositionDataBaseBuilderFILTER_H
+
+// clitk
+#include "clitkFilterWithAnatomicalFeatureDatabaseManagement.h"
+#include "clitkFilterBase.h"
+#include "clitkRelativePositionAnalyzerFilter.h"
+#include "clitkRelativePositionDataBase.h"
+
+// itk
+#include <itkImageToImageFilter.h>
+
+namespace clitk {
+
+ //--------------------------------------------------------------------
+ /*
+ Analyze the relative position of a Target (mask image) according
+ to some Object, in a given Support. Indicate the main important
+ position of this Target according the Object.
+ */
+ //--------------------------------------------------------------------
+
+ template <class ImageType>
+ class RelativePositionDataBaseBuilderFilter:
+ public virtual FilterBase,
+ public clitk::FilterWithAnatomicalFeatureDatabaseManagement,
+ public itk::ImageToImageFilter<ImageType, ImageType>
+ {
+
+ public:
+ /** Standard class typedefs. */
+ typedef itk::ImageToImageFilter<ImageType, ImageType> Superclass;
+ typedef RelativePositionDataBaseBuilderFilter<ImageType> Self;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ /** Method for creation through the object factory. */
+ itkNewMacro(Self);
+
+ /** Run-time type information (and related methods). */
+ itkTypeMacro(RelativePositionDataBaseBuilderFilter, ImageToImageFilter);
+
+ /** Some convenient typedefs. */
+ typedef typename ImageType::ConstPointer ImageConstPointer;
+ typedef typename ImageType::Pointer ImagePointer;
+ typedef typename ImageType::RegionType RegionType;
+ typedef typename ImageType::PixelType PixelType;
+ typedef typename ImageType::SpacingType SpacingType;
+ typedef typename ImageType::SizeType SizeType;
+ typedef typename ImageType::IndexType IndexType;
+ typedef typename ImageType::PointType PointType;
+
+ /** ImageDimension constants */
+ itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension);
+ FILTERBASE_INIT;
+ typedef itk::Image<float, ImageDimension> FloatImageType;
+
+ // Inputs
+ itkGetConstMacro(SupportName, std::string);
+ itkSetMacro(SupportName, std::string);
+
+ itkGetConstMacro(TargetName, std::string);
+ itkSetMacro(TargetName, std::string);
+
+ void AddObjectName(std::string s) { m_ObjectNames.push_back(s); }
+ std::string & GetObjectName(int i) { return m_ObjectNames[i]; }
+ int GetNumberOfObjects() { return m_ObjectNames.size(); }
+
+ // Options
+ itkGetConstMacro(BackgroundValue, PixelType);
+ itkSetMacro(BackgroundValue, PixelType);
+
+ itkGetConstMacro(ForegroundValue, PixelType);
+ itkSetMacro(ForegroundValue, PixelType);
+
+ itkGetConstMacro(NumberOfBins, int);
+ itkSetMacro(NumberOfBins, int);
+
+ itkGetConstMacro(NumberOfAngles, int);
+ itkSetMacro(NumberOfAngles, int);
+
+ itkGetConstMacro(AreaLossTolerance, double);
+ itkSetMacro(AreaLossTolerance, double);
+
+ // For debug
+ void PrintOptions();
+
+ // I dont want to verify inputs information
+ virtual void VerifyInputInformation() { }
+
+ protected:
+ RelativePositionDataBaseBuilderFilter();
+ virtual ~RelativePositionDataBaseBuilderFilter() {}
+
+ PixelType m_BackgroundValue;
+ PixelType m_ForegroundValue;
+
+ ImagePointer m_Support;
+ ImagePointer m_Object;
+ ImagePointer m_Target;
+
+ std::string m_SupportName;
+ std::string m_TargetName;
+ std::vector<std::string> m_ObjectNames;
+
+ int m_NumberOfAngles;
+ std::vector<double> m_ListOfAngles;
+ std::vector<clitk::RelativePositionDirectionType> m_ListOfDirections;
+
+ int m_NumberOfBins;
+ double m_AreaLossTolerance;
+
+ virtual void GenerateOutputInformation();
+ virtual void GenerateInputRequestedRegion();
+ virtual void GenerateData();
+
+ typename FloatImageType::Pointer
+ ComputeFuzzyMap(ImageType * object, ImageType * target, double angle);
+
+ void
+ ComputeOptimalThresholds(FloatImageType * map, ImageType * target, int bins, double tolerance,
+ double & threshold, double & reverseThreshold);
+ private:
+ RelativePositionDataBaseBuilderFilter(const Self&); //purposely not implemented
+ void operator=(const Self&); //purposely not implemented
+
+ }; // end class
+ //--------------------------------------------------------------------
+
+} // end namespace clitk
+//--------------------------------------------------------------------
+
+#include "clitkRelativePositionDataBaseBuilderFilter.txx"
+
+#endif
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+//--------------------------------------------------------------------
+template <class ImageType>
+clitk::RelativePositionDataBaseBuilderFilter<ImageType>::
+RelativePositionDataBaseBuilderFilter():
+ clitk::FilterBase(),
+ clitk::FilterWithAnatomicalFeatureDatabaseManagement(),
+ itk::ImageToImageFilter<ImageType, ImageType>()
+{
+ this->SetNumberOfRequiredInputs(0); // support
+ VerboseFlagOff();
+ SetBackgroundValue(0);
+ SetForegroundValue(1);
+ SetNumberOfBins(100);
+ SetNumberOfAngles(4);
+ SetAreaLossTolerance(0.01);
+ m_ListOfAngles.clear();
+ // SetSupportSize(0);
+ // SetTargetSize(0);
+ // SetSizeWithThreshold(0);
+ // SetSizeWithReverseThreshold(0);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionDataBaseBuilderFilter<ImageType>::
+PrintOptions()
+{
+ DD("TODO");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionDataBaseBuilderFilter<ImageType>::
+GenerateOutputInformation()
+{
+ // ImagePointer input = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ // ImagePointer outputImage = this->GetOutput(0);
+ // outputImage->SetRegions(outputImage->GetLargestPossibleRegion());
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionDataBaseBuilderFilter<ImageType>::
+GenerateInputRequestedRegion()
+{
+ // Call default
+ // itk::ImageToImageFilter<ImageType, ImageType>::GenerateInputRequestedRegion();
+ // // Get input pointers and set requested region to common region
+ // ImagePointer input1 = dynamic_cast<ImageType*>(itk::ProcessObject::GetInput(0));
+ // input1->SetRequestedRegion(input1->GetLargestPossibleRegion());
+}
+//--------------------------------------------------------------------
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::RelativePositionDataBaseBuilderFilter<ImageType>::
+GenerateData()
+{
+ // Load database of anatomical elements
+ static const unsigned int dim = ImageType::ImageDimension;
+ this->LoadAFDB();
+
+ // Get some information
+ std::string patient = this->GetAFDB()->GetTagValue("PatientID");
+
+ // Get input pointers
+ m_Support = this->GetAFDB()->template GetImage <ImageType>(GetSupportName());
+ m_Target = this->GetAFDB()->template GetImage <ImageType>(GetTargetName());
+
+ // Build the list of tested directions
+ m_ListOfAngles.clear();
+ for(uint i=0; i<GetNumberOfAngles(); i++) {
+ double a = i*360.0/GetNumberOfAngles();
+ if (a>180) a = 180-a;
+ m_ListOfAngles.push_back(clitk::deg2rad(a));
+ RelativePositionDirectionType r;
+ r.angle1 = clitk::deg2rad(a);
+ r.angle2 = 0;
+ r.notFlag = false;
+ m_ListOfDirections.push_back(r); // only one direction
+ }
+
+
+ // Perform the RelativePositionAnalyzerFilter for each objects
+ typedef typename clitk::RelativePositionAnalyzerFilter<ImageType> FilterType;
+ for (int i=0; i<GetNumberOfObjects(); i++) {
+ m_Object = this->GetAFDB()->template GetImage <ImageType>(GetObjectName(i));
+
+ for (int j=0; j<m_ListOfDirections.size(); j++) {
+ clitk::RelativePositionDirectionType direction = m_ListOfDirections[j];
+
+ // Create the filter
+ typename FilterType::Pointer filter = FilterType::New();
+ filter->SetInputSupport(m_Support);
+ filter->SetInputTarget(m_Target);
+ filter->SetInputObject(m_Object); // FIXME do AndNot before + only compute supportSize once.
+ filter->SetNumberOfBins(GetNumberOfBins());
+ filter->SetAreaLossTolerance(GetAreaLossTolerance());
+ filter->SetDirection(direction);
+ filter->Update();
+
+ // Results
+ std::ostringstream s;
+ s << patient << " "
+ << GetSupportName() << " "
+ // << GetTargetName() << " " // No need
+ << GetObjectName(i) <<" ";
+ // Normal
+ // if (filter->GetInfo().sizeAfterThreshold != filter->GetInfo().sizeBeforeThreshold) {
+ std::ostringstream os;
+ os << s.str();
+ direction.notFlag = false;
+ direction.Print(os);
+ filter->GetInfo().Print(os);
+ std::cout << os.str() << std::endl;
+ // }
+ // Inverse
+ // if (filter->GetInfoReverse().sizeAfterThreshold != filter->GetInfoReverse().sizeBeforeThreshold) {
+ std::ostringstream oos;
+ oos << s.str();
+ direction.notFlag = true;
+ direction.Print(oos);
+ filter->GetInfoReverse().Print(oos);
+ std::cout << oos.str() << std::endl;
+ // }
+ } // end direction
+
+ } // end object
+}
+//--------------------------------------------------------------------
+
+
namespace clitk {
//--------------------------------------------------------------------
- template<class ImageType>
- void ComputeBBFromImageRegion(const ImageType * image,
- typename ImageType::RegionType region,
- typename itk::BoundingBox<unsigned long,
- ImageType::ImageDimension>::Pointer bb);
-
- //--------------------------------------------------------------------
- template<int Dimension>
- void ComputeBBIntersection(typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbo,
- typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi1,
- typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi2);
-
- //--------------------------------------------------------------------
- template<class ImageType>
- void ComputeRegionFromBB(const ImageType * image,
- const typename itk::BoundingBox<unsigned long,
- ImageType::ImageDimension>::Pointer bb,
- typename ImageType::RegionType & region);
- //--------------------------------------------------------------------
template<class TInternalImageType, class TMaskInternalImageType>
typename TInternalImageType::Pointer
SetBackground(const TInternalImageType * input,
int minimalComponentSize,
LabelizeParameters<typename TImageType::PixelType> * param);
- //--------------------------------------------------------------------
- template<class ImageType>
- typename ImageType::Pointer
- ResizeImageLike(const ImageType * input,
- const itk::ImageBase<ImageType::ImageDimension> * like,
- typename ImageType::PixelType BG);
-
//--------------------------------------------------------------------
template<class MaskImageType>
double spacing=-1,
bool autocropflag=true,
bool singleObjectCCL=true);
+ template<class MaskImageType>
+ typename MaskImageType::Pointer
+ SliceBySliceRelativePosition(const MaskImageType * input,
+ const MaskImageType * object,
+ int direction,
+ double threshold,
+ double angle,
+ bool inverseflag,
+ bool uniqueConnectedComponent=false,
+ double spacing=-1,
+ bool autocropflag=true,
+ bool singleObjectCCL=true);
//--------------------------------------------------------------------
// In a binary image, search for the point belonging to the FG that
namespace clitk {
- //--------------------------------------------------------------------
- template<class ImageType>
- void ComputeBBFromImageRegion(const ImageType * image,
- typename ImageType::RegionType region,
- typename itk::BoundingBox<unsigned long,
- ImageType::ImageDimension>::Pointer bb) {
- typedef typename ImageType::IndexType IndexType;
- IndexType firstIndex;
- IndexType lastIndex;
- for(unsigned int i=0; i<image->GetImageDimension(); i++) {
- firstIndex[i] = region.GetIndex()[i];
- lastIndex[i] = firstIndex[i]+region.GetSize()[i];
- }
-
- typedef itk::BoundingBox<unsigned long,
- ImageType::ImageDimension> BBType;
- typedef typename BBType::PointType PointType;
- PointType lastPoint;
- PointType firstPoint;
- image->TransformIndexToPhysicalPoint(firstIndex, firstPoint);
- image->TransformIndexToPhysicalPoint(lastIndex, lastPoint);
-
- bb->SetMaximum(lastPoint);
- bb->SetMinimum(firstPoint);
- }
- //--------------------------------------------------------------------
-
-
- //--------------------------------------------------------------------
- template<int Dimension>
- void ComputeBBIntersection(typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbo,
- typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi1,
- typename itk::BoundingBox<unsigned long, Dimension>::Pointer bbi2) {
-
- typedef itk::BoundingBox<unsigned long, Dimension> BBType;
- typedef typename BBType::PointType PointType;
- PointType lastPoint;
- PointType firstPoint;
-
- for(unsigned int i=0; i<Dimension; i++) {
- firstPoint[i] = std::max(bbi1->GetMinimum()[i],
- bbi2->GetMinimum()[i]);
- lastPoint[i] = std::min(bbi1->GetMaximum()[i],
- bbi2->GetMaximum()[i]);
- }
-
- bbo->SetMaximum(lastPoint);
- bbo->SetMinimum(firstPoint);
- }
- //--------------------------------------------------------------------
-
-
- //--------------------------------------------------------------------
- template<class ImageType>
- void ComputeRegionFromBB(const ImageType * image,
- const typename itk::BoundingBox<unsigned long,
- ImageType::ImageDimension>::Pointer bb,
- typename ImageType::RegionType & region) {
- // Types
- typedef typename ImageType::IndexType IndexType;
- typedef typename ImageType::PointType PointType;
- typedef typename ImageType::RegionType RegionType;
- typedef typename ImageType::SizeType SizeType;
-
- // Region starting point
- IndexType regionStart;
- PointType start = bb->GetMinimum();
- image->TransformPhysicalPointToIndex(start, regionStart);
-
- // Region size
- SizeType regionSize;
- PointType maxs = bb->GetMaximum();
- PointType mins = bb->GetMinimum();
- for(unsigned int i=0; i<ImageType::ImageDimension; i++) {
- regionSize[i] = lrint((maxs[i] - mins[i])/image->GetSpacing()[i]);
- }
-
- // Create region
- region.SetIndex(regionStart);
- region.SetSize(regionSize);
- }
- //--------------------------------------------------------------------
-
//--------------------------------------------------------------------
template<class ImageType, class TMaskImageType>
typename ImageType::Pointer
//--------------------------------------------------------------------
- template<class ImageType>
- typename ImageType::Pointer
- ResizeImageLike(const ImageType * input,
- const itk::ImageBase<ImageType::ImageDimension> * like,
- typename ImageType::PixelType backgroundValue)
+ template<class MaskImageType>
+ typename MaskImageType::Pointer
+ SliceBySliceRelativePosition(const MaskImageType * input,
+ const MaskImageType * object,
+ int direction,
+ double threshold,
+ std::string orientation,
+ bool uniqueConnectedComponent,
+ double spacing,
+ bool autocropFlag,
+ bool singleObjectCCL)
{
- typedef CropLikeImageFilter<ImageType> CropFilterType;
- typename CropFilterType::Pointer cropFilter = CropFilterType::New();
- cropFilter->SetInput(input);
- cropFilter->SetCropLikeImage(like);
- cropFilter->SetBackgroundValue(backgroundValue);
- cropFilter->Update();
- return cropFilter->GetOutput();
+ typedef clitk::SliceBySliceRelativePositionFilter<MaskImageType> SliceRelPosFilterType;
+ typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New();
+ sliceRelPosFilter->VerboseStepFlagOff();
+ sliceRelPosFilter->WriteStepFlagOff();
+ sliceRelPosFilter->SetInput(input);
+ sliceRelPosFilter->SetInputObject(object);
+ sliceRelPosFilter->SetDirection(direction);
+ sliceRelPosFilter->SetFuzzyThreshold(threshold);
+ sliceRelPosFilter->AddOrientationTypeString(orientation);
+ sliceRelPosFilter->SetIntermediateSpacingFlag((spacing != -1));
+ sliceRelPosFilter->SetIntermediateSpacing(spacing);
+ sliceRelPosFilter->SetUniqueConnectedComponentBySliceFlag(uniqueConnectedComponent);
+ sliceRelPosFilter->ObjectCCLSelectionFlagOff();
+ sliceRelPosFilter->SetUseTheLargestObjectCCLFlag(singleObjectCCL);
+ // sliceRelPosFilter->SetInverseOrientationFlag(inverseflag);
+ sliceRelPosFilter->SetAutoCropFlag(autocropFlag);
+ sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn();
+ sliceRelPosFilter->Update();
+ return sliceRelPosFilter->GetOutput();
}
//--------------------------------------------------------------------
const MaskImageType * object,
int direction,
double threshold,
- std::string orientation,
+ double angle,
+ bool inverseflag,
bool uniqueConnectedComponent,
double spacing,
bool autocropFlag,
bool singleObjectCCL)
{
- typedef SliceBySliceRelativePositionFilter<MaskImageType> SliceRelPosFilterType;
+ typedef clitk::SliceBySliceRelativePositionFilter<MaskImageType> SliceRelPosFilterType;
typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New();
sliceRelPosFilter->VerboseStepFlagOff();
sliceRelPosFilter->WriteStepFlagOff();
sliceRelPosFilter->SetInputObject(object);
sliceRelPosFilter->SetDirection(direction);
sliceRelPosFilter->SetFuzzyThreshold(threshold);
- sliceRelPosFilter->AddOrientationTypeString(orientation);
+ // sliceRelPosFilter->AddOrientationTypeString(orientation);
+ sliceRelPosFilter->AddAnglesInRad(angle, 0.0);
sliceRelPosFilter->SetIntermediateSpacingFlag((spacing != -1));
sliceRelPosFilter->SetIntermediateSpacing(spacing);
sliceRelPosFilter->SetUniqueConnectedComponentBySliceFlag(uniqueConnectedComponent);
sliceRelPosFilter->ObjectCCLSelectionFlagOff();
sliceRelPosFilter->SetUseTheLargestObjectCCLFlag(singleObjectCCL);
- // sliceRelPosFilter->SetInverseOrientationFlag(inverseflag);
+ sliceRelPosFilter->SetInverseOrientationFlag(inverseflag);
sliceRelPosFilter->SetAutoCropFlag(autocropFlag);
sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn();
sliceRelPosFilter->Update();
typename ImageType::PointType p;
image->TransformIndexToPhysicalPoint(image->GetLargestPossibleRegion().GetIndex()+
image->GetLargestPossibleRegion().GetSize(), p);
- return CropImageAlongOneAxis<ImageType>(image, dim, max, p[dim], autoCrop, BG);
+ // Add GetSpacing because remove Lower or equal than
+ // DD(max);
+ // DD(p);
+ // DD(max+image->GetSpacing()[dim]);
+ return CropImageAlongOneAxis<ImageType>(image, dim, max+image->GetSpacing()[dim], p[dim], autoCrop, BG);
}
//--------------------------------------------------------------------
/** ImageDimension constants */
itkStaticConstMacro(ImageDimension, unsigned int, ImageType::ImageDimension);
typedef itk::Image<float, ImageDimension> FloatImageType;
+ typedef itk::Image<float, ImageDimension-1> FloatSliceType;
/** Some convenient typedefs. */
typedef typename ImageType::ConstPointer ImageConstPointer;
======================================================================-====*/
// clitk
+#include "clitkCropLikeImageFilter.h"
#include "clitkSegmentationUtils.h"
#include "clitkExtractSliceFilter.h"
#include "clitkResampleImageWithOptionsFilter.h"
{
os << "Slice direction = " << this->GetDirection() << std::endl
<< "BG value = " << this->GetBackgroundValue() << std::endl;
- for(int i=0; i<this->GetNumberOfAngles(); i++)
+ for(int i=0; i<this->GetNumberOfAngles(); i++) {
os << "Orientation = " << this->GetOrientationTypeString()[i] << std::endl;
+ os << "Angles = " << clitk::rad2deg(this->GetAngle1InRad(i))
+ << " " << clitk::rad2deg(this->GetAngle2InRad(i)) << std::endl;
+ }
os << "InverseOrientationFlag = " << this->GetInverseOrientationFlag() << std::endl
<< "SpacingFlag = " << this->GetIntermediateSpacingFlag() << std::endl
<< "Spacing = " << this->GetIntermediateSpacing() << std::endl
extractSliceFilter->GetOutputSlices(mObjectSlices);
this->template StopCurrentStep<SliceType>(mObjectSlices[0]);
+ //--------------------------------------------------------------------
+ // Prepare fuzzy slices (if needed)
+ std::vector<typename FloatSliceType::Pointer> mFuzzyMapSlices;
+ mFuzzyMapSlices.resize(mInputSlices.size());
+
//--------------------------------------------------------------------
// Perform slice by slice relative position
this->StartNewStep("Perform slice by slice relative position");
// Count the number of CCL (allow to ignore empty slice)
int nb=0;
mObjectSlices[i] = LabelizeAndCountNumberOfObjects<SliceType>(mObjectSlices[i], 0, true, 1, nb);
- if ((!GetIgnoreEmptySliceObjectFlag()) || (nb!=0)) {
- // Select or not a single CCL ?
- if (GetUseTheLargestObjectCCLFlag()) {
- mObjectSlices[i] = KeepLabels<SliceType>(mObjectSlices[i], 0, 1, 1, 1, true);
- }
+ // If no object and empty slices :
+ if ((nb==0) && (this->GetFuzzyMapOnlyFlag())) {
+ typename FloatSliceType::Pointer one = FloatSliceType::New();
+ one->CopyInformation(mObjectSlices[0]);
+ one->SetRegions(mObjectSlices[0]->GetLargestPossibleRegion());
+ one->Allocate();
+ one->FillBuffer(2.0);
+ mFuzzyMapSlices[i] = one;
+ }
+ else {
+ if ((!GetIgnoreEmptySliceObjectFlag()) || (nb!=0)) {
- // Select a single according to a position if more than one CCL
- if (GetObjectCCLSelectionFlag()) {
- // if several CCL, choose the most extrema according a direction,
- // if not -> should we consider this slice ?
- if (nb<2) {
- if (GetObjectCCLSelectionIgnoreSingleCCLFlag()) {
- mObjectSlices[i] = SetBackground<SliceType, SliceType>(mObjectSlices[i], mObjectSlices[i],
- 1, this->GetBackgroundValue(),
- true);
- }
+ // Select or not a single CCL ?
+ if (GetUseTheLargestObjectCCLFlag()) {
+ mObjectSlices[i] = KeepLabels<SliceType>(mObjectSlices[i], 0, 1, 1, 1, true);
}
- int dim = GetObjectCCLSelectionDimension();
- int direction = GetObjectCCLSelectionDirection();
- std::vector<typename SliceType::PointType> centroids;
- ComputeCentroids<SliceType>(mObjectSlices[i], this->GetBackgroundValue(), centroids);
- uint index=1;
- for(uint j=1; j<centroids.size(); j++) {
- if (direction == 1) {
- if (centroids[j][dim] > centroids[index][dim]) index = j;
+
+ // Select a single according to a position if more than one CCL
+ if (GetObjectCCLSelectionFlag()) {
+ // if several CCL, choose the most extrema according a direction,
+ // if not -> should we consider this slice ?
+ if (nb<2) {
+ if (GetObjectCCLSelectionIgnoreSingleCCLFlag()) {
+ mObjectSlices[i] = SetBackground<SliceType, SliceType>(mObjectSlices[i], mObjectSlices[i],
+ 1, this->GetBackgroundValue(),
+ true);
+ }
+ }
+ int dim = GetObjectCCLSelectionDimension();
+ int direction = GetObjectCCLSelectionDirection();
+ std::vector<typename SliceType::PointType> centroids;
+ ComputeCentroids<SliceType>(mObjectSlices[i], this->GetBackgroundValue(), centroids);
+ uint index=1;
+ for(uint j=1; j<centroids.size(); j++) {
+ if (direction == 1) {
+ if (centroids[j][dim] > centroids[index][dim]) index = j;
+ }
+ else {
+ if (centroids[j][dim] < centroids[index][dim]) index = j;
+ }
}
- else {
- if (centroids[j][dim] < centroids[index][dim]) index = j;
+ for(uint v=1; v<centroids.size(); v++) {
+ if (v != index) {
+ mObjectSlices[i] = SetBackground<SliceType, SliceType>(mObjectSlices[i], mObjectSlices[i],
+ (char)v, this->GetBackgroundValue(),
+ true);
+ }
}
+ } // end GetbjectCCLSelectionFlag = true
+
+ // Relative position
+ typedef clitk::AddRelativePositionConstraintToLabelImageFilter<SliceType> RelPosFilterType;
+ typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
+
+ relPosFilter->VerboseStepFlagOff();
+ relPosFilter->WriteStepFlagOff();
+ relPosFilter->SetBackgroundValue(this->GetBackgroundValue());
+ relPosFilter->SetInput(mInputSlices[i]);
+ relPosFilter->SetInputObject(mObjectSlices[i]);
+ relPosFilter->SetRemoveObjectFlag(this->GetRemoveObjectFlag());
+ // This flag (InverseOrientation) *must* be set before
+ // AddOrientation because AddOrientation can change it.
+ relPosFilter->SetInverseOrientationFlag(this->GetInverseOrientationFlag());
+ for(int j=0; j<this->GetNumberOfAngles(); j++) {
+ // relPosFilter->AddOrientationTypeString(this->GetOrientationTypeString(j));
+ relPosFilter->AddAnglesInRad(this->GetAngle1InRad(j), this->GetAngle2InRad(j));
+ // DD(this->GetOrientationTypeString(j));
}
- for(uint v=1; v<centroids.size(); v++) {
- if (v != index) {
- mObjectSlices[i] = SetBackground<SliceType, SliceType>(mObjectSlices[i], mObjectSlices[i],
- (char)v, this->GetBackgroundValue(),
- true);
+ // DD(this->GetInverseOrientationFlag());
+ //relPosFilter->SetOrientationType(this->GetOrientationType());
+ relPosFilter->SetIntermediateSpacing(this->GetIntermediateSpacing());
+ relPosFilter->SetIntermediateSpacingFlag(this->GetIntermediateSpacingFlag());
+ relPosFilter->SetFuzzyThreshold(this->GetFuzzyThreshold());
+ relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop
+ relPosFilter->SetCombineWithOrFlag(this->GetCombineWithOrFlag());
+
+ // should we stop after fuzzy map ?
+ relPosFilter->SetFuzzyMapOnlyFlag(this->GetFuzzyMapOnlyFlag());
+
+ // Go !
+ relPosFilter->Update();
+
+ // If we stop after the fuzzy map, store the fuzzy slices
+ if (this->GetFuzzyMapOnlyFlag()) {
+ mFuzzyMapSlices[i] = relPosFilter->GetFuzzyMap();
+ // writeImage<FloatSliceType>(mFuzzyMapSlices[i], "slice_"+toString(i)+".mha");
+ }
+ else {
+ mInputSlices[i] = relPosFilter->GetOutput();
+ // Select main CC if needed
+ if (GetUniqueConnectedComponentBySliceFlag()) {
+ mInputSlices[i] = Labelize<SliceType>(mInputSlices[i], 0, true, 1);
+ mInputSlices[i] = KeepLabels<SliceType>(mInputSlices[i], 0, 1, 1, 1, true);
}
+
}
- } // end GetbjectCCLSelectionFlag = true
-
- // Relative position
- typedef clitk::AddRelativePositionConstraintToLabelImageFilter<SliceType> RelPosFilterType;
- typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
-
- relPosFilter->VerboseStepFlagOff();
- relPosFilter->WriteStepFlagOff();
- relPosFilter->SetBackgroundValue(this->GetBackgroundValue());
- relPosFilter->SetInput(mInputSlices[i]);
- relPosFilter->SetInputObject(mObjectSlices[i]);
- relPosFilter->SetRemoveObjectFlag(this->GetRemoveObjectFlag());
- // This flag (InverseOrientation) *must* be set before
- // AddOrientation because AddOrientation can change it.
- relPosFilter->SetInverseOrientationFlag(this->GetInverseOrientationFlag());
- for(int j=0; j<this->GetNumberOfAngles(); j++) {
- relPosFilter->AddOrientationTypeString(this->GetOrientationTypeString(j));
- //DD(this->GetOrientationTypeString(j));
- }
- //DD(this->GetInverseOrientationFlag());
- //relPosFilter->SetOrientationType(this->GetOrientationType());
- relPosFilter->SetIntermediateSpacing(this->GetIntermediateSpacing());
- relPosFilter->SetIntermediateSpacingFlag(this->GetIntermediateSpacingFlag());
- relPosFilter->SetFuzzyThreshold(this->GetFuzzyThreshold());
- relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop
- relPosFilter->SetCombineWithOrFlag(this->GetCombineWithOrFlag());
- relPosFilter->Update();
- mInputSlices[i] = relPosFilter->GetOutput();
-
- // Select main CC if needed
- if (GetUniqueConnectedComponentBySliceFlag()) {
- mInputSlices[i] = Labelize<SliceType>(mInputSlices[i], 0, true, 1);
- mInputSlices[i] = KeepLabels<SliceType>(mInputSlices[i], 0, 1, 1, 1, true);
- }
+ }
/*
// Select unique CC according to the most in a given direction
if (GetUniqueConnectedComponentBySliceAccordingToADirection()) {
- int nb;
- mInputSlices[i] = LabelizeAndCountNumberOfObjects<SliceType>(mInputSlices[i], 0, true, 1, nb);
- std::vector<typename ImageType::PointType> & centroids;
- ComputeCentroids
- }
+ int nb;
+ mInputSlices[i] = LabelizeAndCountNumberOfObjects<SliceType>(mInputSlices[i], 0, true, 1, nb);
+ std::vector<typename ImageType::PointType> & centroids;
+ ComputeCentroids
+ }
*/
}
}
+ // Join the fuzzy map if needed
+ if (this->GetFuzzyMapOnlyFlag()) {
+ this->m_FuzzyMap = clitk::JoinSlices<FloatImageType>(mFuzzyMapSlices, input, GetDirection());
+ this->template StopCurrentStep<FloatImageType>(this->m_FuzzyMap);
+ return;
+ }
+
// Join the slices
m_working_input = clitk::JoinSlices<ImageType>(mInputSlices, input, GetDirection());
this->template StopCurrentStep<ImageType>(m_working_input);
//--------------------------------------------------------------------
// Final Step -> set output
//this->SetNthOutput(0, m_working_input);
+ if (this->GetFuzzyMapOnlyFlag()) return; // no output in this case
this->GraftOutput(m_working_input);
return;
}
--- /dev/null
+/*=========================================================================
+
+ Program: Insight Segmentation & Registration Toolkit
+ Module: $RCSfile: itkLabelOverlapMeasuresImageFilter.h,v $
+ Language: C++
+ Date: $Date: $
+ Version: $Revision: $
+
+ Copyright (c) Insight Software Consortium. All rights reserved.
+ See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the above copyright notices for more information.
+
+=========================================================================*/
+#ifndef __itkLabelOverlapMeasuresImageFilter_h
+#define __itkLabelOverlapMeasuresImageFilter_h
+
+#include "itkImageToImageFilter.h"
+#include "itkFastMutexLock.h"
+#include "itkNumericTraits.h"
+
+//#include "itk_hash_map.h"
+#include "itksys/hash_map.hxx"
+
+namespace itk {
+
+/** \class LabelOverlapMeasuresImageFilter
+ * \brief Computes overlap measures between the set same set of labels of
+ * pixels of two images. Background is assumed to be 0.
+ *
+ * \sa LabelOverlapMeasuresImageFilter
+ *
+ * \ingroup MultiThreaded
+ */
+template<class TLabelImage>
+class ITK_EXPORT LabelOverlapMeasuresImageFilter :
+ public ImageToImageFilter<TLabelImage, TLabelImage>
+{
+public:
+ /** Standard Self typedef */
+ typedef LabelOverlapMeasuresImageFilter Self;
+ typedef ImageToImageFilter<TLabelImage,TLabelImage> Superclass;
+ typedef SmartPointer<Self> Pointer;
+ typedef SmartPointer<const Self> ConstPointer;
+
+ /** Method for creation through the object factory. */
+ itkNewMacro( Self );
+
+ /** Runtime information support. */
+ itkTypeMacro( LabelOverlapMeasuresImageFilter, ImageToImageFilter );
+
+ virtual void VerifyInputInformation() { }
+
+
+ /** Image related typedefs. */
+ typedef TLabelImage LabelImageType;
+ typedef typename TLabelImage::Pointer LabelImagePointer;
+ typedef typename TLabelImage::ConstPointer LabelImageConstPointer;
+
+ typedef typename TLabelImage::RegionType RegionType;
+ typedef typename TLabelImage::SizeType SizeType;
+ typedef typename TLabelImage::IndexType IndexType;
+
+ typedef typename TLabelImage::PixelType LabelType;
+
+ /** Type to use form computations. */
+ typedef typename NumericTraits<LabelType>::RealType RealType;
+
+ /** \class LabelLabelOverlapMeasuress
+ * \brief Metrics stored per label */
+ class LabelSetMeasures
+ {
+ public:
+ // default constructor
+ LabelSetMeasures()
+ {
+ m_Source = 0;
+ m_Target = 0;
+ m_Union = 0;
+ m_Intersection = 0;
+ m_SourceComplement = 0;
+ m_TargetComplement = 0;
+ }
+
+ // added for completeness
+ LabelSetMeasures& operator=( const LabelSetMeasures& l )
+ {
+ m_Source = l.m_Source;
+ m_Target = l.m_Target;
+ m_Union = l.m_Union;
+ m_Intersection = l.m_Intersection;
+ m_SourceComplement = l.m_SourceComplement;
+ m_TargetComplement = l.m_TargetComplement;
+ }
+
+ unsigned long m_Source;
+ unsigned long m_Target;
+ unsigned long m_Union;
+ unsigned long m_Intersection;
+ unsigned long m_SourceComplement;
+ unsigned long m_TargetComplement;
+ };
+
+ /** Type of the map used to store data per label */
+ typedef itksys::hash_map<LabelType, LabelSetMeasures> MapType;
+ typedef typename MapType::iterator MapIterator;
+ typedef typename MapType::const_iterator MapConstIterator;
+
+ /** Image related typedefs. */
+ itkStaticConstMacro( ImageDimension, unsigned int,
+ TLabelImage::ImageDimension );
+
+ /** Set the source image. */
+ void SetSourceImage( const LabelImageType * image )
+ { this->SetNthInput( 0, const_cast<LabelImageType *>( image ) ); }
+
+ /** Set the target image. */
+ void SetTargetImage( const LabelImageType * image )
+ { this->SetNthInput( 1, const_cast<LabelImageType *>( image ) ); }
+
+ /** Get the source image. */
+ const LabelImageType * GetSourceImage( void )
+ { return this->GetInput( 0 ); }
+
+ /** Get the target image. */
+ const LabelImageType * GetTargetImage( void )
+ { return this->GetInput( 1 ); }
+
+ /** Get the label set measures */
+ MapType GetLabelSetMeasures()
+ { return this->m_LabelSetMeasures; }
+
+ /**
+ * tric overlap measures
+ */
+ /** measures over all labels */
+ RealType GetTotalOverlap();
+ RealType GetUnionOverlap();
+ RealType GetMeanOverlap();
+ RealType GetVolumeSimilarity();
+ RealType GetFalseNegativeError();
+ RealType GetFalsePositiveError();
+ /** measures over individual labels */
+ RealType GetTargetOverlap( LabelType );
+ RealType GetUnionOverlap( LabelType );
+ RealType GetMeanOverlap( LabelType );
+ RealType GetVolumeSimilarity( LabelType );
+ RealType GetFalseNegativeError( LabelType );
+ RealType GetFalsePositiveError( LabelType );
+ /** alternative names */
+ RealType GetJaccardCoefficient()
+ { return this->GetUnionOverlap(); }
+ RealType GetJaccardCoefficient( LabelType label )
+ { return this->GetUnionOverlap( label ); }
+ RealType GetDiceCoefficient()
+ { return this->GetMeanOverlap(); }
+ RealType GetDiceCoefficient( LabelType label )
+ { return this->GetMeanOverlap( label ); }
+
+
+#ifdef ITK_USE_CONCEPT_CHECKING
+ /** Begin concept checking */
+ itkConceptMacro( Input1HasNumericTraitsCheck,
+ ( Concept::HasNumericTraits<LabelType> ) );
+ /** End concept checking */
+#endif
+
+protected:
+ LabelOverlapMeasuresImageFilter();
+ ~LabelOverlapMeasuresImageFilter(){};
+ void PrintSelf( std::ostream& os, Indent indent ) const;
+
+ /**
+ * Pass the input through unmodified. Do this by setting the output to the
+ * source this by setting the output to the source image in the
+ * AllocateOutputs() method.
+ */
+ void AllocateOutputs();
+
+ void BeforeThreadedGenerateData();
+
+ void AfterThreadedGenerateData();
+
+ /** Multi-thread version GenerateData. */
+ void ThreadedGenerateData( const RegionType&, int );
+
+ // Override since the filter needs all the data for the algorithm
+ void GenerateInputRequestedRegion();
+
+ // Override since the filter produces all of its output
+ void EnlargeOutputRequestedRegion( DataObject *data );
+
+private:
+ LabelOverlapMeasuresImageFilter( const Self& ); //purposely not implemented
+ void operator=( const Self& ); //purposely not implemented
+
+ std::vector<MapType> m_LabelSetMeasuresPerThread;
+ MapType m_LabelSetMeasures;
+
+ SimpleFastMutexLock m_Mutex;
+
+}; // end of class
+
+} // end namespace itk
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include "itkLabelOverlapMeasuresImageFilter.txx"
+#endif
+
+#endif
--- /dev/null
+/*=========================================================================
+
+ Program: Insight Segmentation & Registration Toolkit
+ Module: $RCSfile: itkLabelOverlapMeasuresImageFilter.txx,v $
+ Language: C++
+ Date: $Date: $
+ Version: $Revision: $
+
+ Copyright (c) Insight Software Consortium. All rights reserved.
+ See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the above copyright notices for more information.
+
+=========================================================================*/
+#ifndef _itkLabelOverlapMeasuresImageFilter_txx
+#define _itkLabelOverlapMeasuresImageFilter_txx
+
+#include "itkLabelOverlapMeasuresImageFilter.h"
+
+#include "itkImageRegionConstIterator.h"
+#include "itkProgressReporter.h"
+
+namespace itk {
+
+#if defined(__GNUC__) && (__GNUC__ <= 2) //NOTE: This class needs a mutex for gnu 2.95
+/** Used for mutex locking */
+#define LOCK_HASHMAP this->m_Mutex.Lock()
+#define UNLOCK_HASHMAP this->m_Mutex.Unlock()
+#else
+#define LOCK_HASHMAP
+#define UNLOCK_HASHMAP
+#endif
+
+template<class TLabelImage>
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::LabelOverlapMeasuresImageFilter()
+{
+ // this filter requires two input images
+ this->SetNumberOfRequiredInputs( 2 );
+}
+
+template<class TLabelImage>
+void
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GenerateInputRequestedRegion()
+{
+ Superclass::GenerateInputRequestedRegion();
+ if( this->GetSourceImage() )
+ {
+ LabelImagePointer source = const_cast
+ <LabelImageType *>( this->GetSourceImage() );
+ source->SetRequestedRegionToLargestPossibleRegion();
+ }
+ if( this->GetTargetImage() )
+ {
+ LabelImagePointer target = const_cast
+ <LabelImageType *>( this->GetTargetImage() );
+ target->SetRequestedRegionToLargestPossibleRegion();
+ }
+}
+
+template<class TLabelImage>
+void
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::EnlargeOutputRequestedRegion( DataObject *data )
+{
+ Superclass::EnlargeOutputRequestedRegion( data );
+ data->SetRequestedRegionToLargestPossibleRegion();
+}
+
+
+template<class TLabelImage>
+void
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::AllocateOutputs()
+{
+ // Pass the source through as the output
+ LabelImagePointer image =
+ const_cast<TLabelImage *>( this->GetSourceImage() );
+ this->SetNthOutput( 0, image );
+
+ // Nothing that needs to be allocated for the remaining outputs
+}
+
+template<class TLabelImage>
+void
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::BeforeThreadedGenerateData()
+{
+ int numberOfThreads = this->GetNumberOfThreads();
+
+ // Resize the thread temporaries
+ this->m_LabelSetMeasuresPerThread.resize( numberOfThreads );
+
+ // Initialize the temporaries
+ for( int n = 0; n < numberOfThreads; n++ )
+ {
+ this->m_LabelSetMeasuresPerThread[n].clear();
+ }
+
+ // Initialize the final map
+ this->m_LabelSetMeasures.clear();
+}
+
+template<class TLabelImage>
+void
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::AfterThreadedGenerateData()
+{
+ // Run through the map for each thread and accumulate the set measures.
+ for( int n = 0; n < this->GetNumberOfThreads(); n++ )
+ {
+ // iterate over the map for this thread
+ for( MapConstIterator threadIt = this->m_LabelSetMeasuresPerThread[n].begin();
+ threadIt != this->m_LabelSetMeasuresPerThread[n].end();
+ ++threadIt )
+ {
+ // does this label exist in the cumulative stucture yet?
+ MapIterator mapIt = this->m_LabelSetMeasures.find( ( *threadIt ).first );
+ if( mapIt == this->m_LabelSetMeasures.end() )
+ {
+ // create a new entry
+ typedef typename MapType::value_type MapValueType;
+ mapIt = this->m_LabelSetMeasures.insert( MapValueType(
+ (*threadIt).first, LabelSetMeasures() ) ).first;
+ }
+
+ // accumulate the information from this thread
+ (*mapIt).second.m_Source += (*threadIt).second.m_Source;
+ (*mapIt).second.m_Target += (*threadIt).second.m_Target;
+ (*mapIt).second.m_Union += (*threadIt).second.m_Union;
+ (*mapIt).second.m_Intersection +=
+ (*threadIt).second.m_Intersection;
+ (*mapIt).second.m_SourceComplement +=
+ (*threadIt).second.m_SourceComplement;
+ (*mapIt).second.m_TargetComplement +=
+ (*threadIt).second.m_TargetComplement;
+ } // end of thread map iterator loop
+ } // end of thread loop
+}
+
+template<class TLabelImage>
+void
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::ThreadedGenerateData( const RegionType& outputRegionForThread,
+ int threadId )
+{
+ ImageRegionConstIterator<LabelImageType> ItS( this->GetSourceImage(),
+ outputRegionForThread );
+ ImageRegionConstIterator<LabelImageType> ItT( this->GetTargetImage(),
+ outputRegionForThread );
+
+ // support progress methods/callbacks
+ ProgressReporter progress( this, threadId,
+ 2*outputRegionForThread.GetNumberOfPixels() );
+
+ for( ItS.GoToBegin(), ItT.GoToBegin(); !ItS.IsAtEnd(); ++ItS, ++ItT )
+ {
+ LabelType sourceLabel = ItS.Get();
+ LabelType targetLabel = ItT.Get();
+
+ // is the label already in this thread?
+ MapIterator mapItS =
+ this->m_LabelSetMeasuresPerThread[threadId].find( sourceLabel );
+ MapIterator mapItT =
+ this->m_LabelSetMeasuresPerThread[threadId].find( targetLabel );
+
+ if( mapItS == this->m_LabelSetMeasuresPerThread[threadId].end() )
+ {
+ // create a new label set measures object
+ typedef typename MapType::value_type MapValueType;
+ LOCK_HASHMAP;
+ mapItS = this->m_LabelSetMeasuresPerThread[threadId].insert(
+ MapValueType( sourceLabel, LabelSetMeasures() ) ).first;
+ UNLOCK_HASHMAP;
+ }
+
+ if( mapItT == this->m_LabelSetMeasuresPerThread[threadId].end() )
+ {
+ // create a new label set measures object
+ typedef typename MapType::value_type MapValueType;
+ LOCK_HASHMAP;
+ mapItT = this->m_LabelSetMeasuresPerThread[threadId].insert(
+ MapValueType( targetLabel, LabelSetMeasures() ) ).first;
+ UNLOCK_HASHMAP;
+ }
+
+ (*mapItS).second.m_Source++;
+ (*mapItT).second.m_Target++;
+
+ if( sourceLabel == targetLabel )
+ {
+ (*mapItS).second.m_Intersection++;
+ (*mapItS).second.m_Union++;
+ }
+ else
+ {
+ (*mapItS).second.m_Union++;
+ (*mapItT).second.m_Union++;
+
+ (*mapItS).second.m_SourceComplement++;
+ (*mapItT).second.m_TargetComplement++;
+ }
+
+ progress.CompletedPixel();
+ }
+}
+
+/**
+ * measures
+ */
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetTotalOverlap()
+{
+ RealType numerator = 0.0;
+ RealType denominator = 0.0;
+ for( MapIterator mapIt = this->m_LabelSetMeasures.begin();
+ mapIt != this->m_LabelSetMeasures.end(); ++mapIt )
+ {
+ // Do not include the background in the final value.
+ if( (*mapIt).first == NumericTraits<LabelType>::Zero )
+ {
+ continue;
+ }
+ numerator += static_cast<RealType>( (*mapIt).second.m_Intersection );
+ denominator += static_cast<RealType>( (*mapIt).second.m_Target );
+ }
+ return ( numerator / denominator );
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetTargetOverlap( LabelType label )
+{
+ MapIterator mapIt = this->m_LabelSetMeasures.find( label );
+ if( mapIt == this->m_LabelSetMeasures.end() )
+ {
+ itkWarningMacro( "Label " << label << " not found." );
+ return 0.0;
+ }
+ RealType value =
+ static_cast<RealType>( (*mapIt).second.m_Intersection ) /
+ static_cast<RealType>( (*mapIt).second.m_Target );
+ return value;
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetUnionOverlap()
+{
+ RealType numerator = 0.0;
+ RealType denominator = 0.0;
+ for( MapIterator mapIt = this->m_LabelSetMeasures.begin();
+ mapIt != this->m_LabelSetMeasures.end(); ++mapIt )
+ {
+ // Do not include the background in the final value.
+ if( (*mapIt).first == NumericTraits<LabelType>::Zero )
+ {
+ continue;
+ }
+ numerator += static_cast<RealType>( (*mapIt).second.m_Intersection );
+ denominator += static_cast<RealType>( (*mapIt).second.m_Union );
+ }
+ return ( numerator / denominator );
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetUnionOverlap( LabelType label )
+{
+ MapIterator mapIt = this->m_LabelSetMeasures.find( label );
+ if( mapIt == this->m_LabelSetMeasures.end() )
+ {
+ itkWarningMacro( "Label " << label << " not found." );
+ return 0.0;
+ }
+ RealType value =
+ static_cast<RealType>( (*mapIt).second.m_Intersection ) /
+ static_cast<RealType>( (*mapIt).second.m_Union );
+ return value;
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetMeanOverlap()
+{
+ RealType uo = this->GetUnionOverlap();
+ return ( 2.0 * uo / ( 1.0 + uo ) );
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetMeanOverlap( LabelType label )
+{
+ RealType uo = this->GetUnionOverlap( label );
+ return ( 2.0 * uo / ( 1.0 + uo ) );
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetVolumeSimilarity()
+{
+ RealType numerator = 0.0;
+ RealType denominator = 0.0;
+ for( MapIterator mapIt = this->m_LabelSetMeasures.begin();
+ mapIt != this->m_LabelSetMeasures.end(); ++mapIt )
+ {
+ // Do not include the background in the final value.
+ if( (*mapIt).first == NumericTraits<LabelType>::Zero )
+ {
+ continue;
+ }
+ numerator += ( ( static_cast<RealType>( (*mapIt).second.m_Source ) -
+ static_cast<RealType>( (*mapIt).second.m_Target ) ) );
+ denominator += ( ( static_cast<RealType>( (*mapIt).second.m_Source ) +
+ static_cast<RealType>( (*mapIt).second.m_Target ) ) );
+ }
+ return ( 2.0 * numerator / denominator );
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetVolumeSimilarity( LabelType label )
+{
+ MapIterator mapIt = this->m_LabelSetMeasures.find( label );
+ if( mapIt == this->m_LabelSetMeasures.end() )
+ {
+ itkWarningMacro( "Label " << label << " not found." );
+ return 0.0;
+ }
+ RealType value = 2.0 *
+ ( static_cast<RealType>( (*mapIt).second.m_Source ) -
+ static_cast<RealType>( (*mapIt).second.m_Target ) ) /
+ ( static_cast<RealType>( (*mapIt).second.m_Source ) +
+ static_cast<RealType>( (*mapIt).second.m_Target ) );
+ return value;
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetFalseNegativeError()
+{
+ RealType numerator = 0.0;
+ RealType denominator = 0.0;
+ for( MapIterator mapIt = this->m_LabelSetMeasures.begin();
+ mapIt != this->m_LabelSetMeasures.end(); ++mapIt )
+ {
+ // Do not include the background in the final value.
+ if( (*mapIt).first == NumericTraits<LabelType>::Zero )
+ {
+ continue;
+ }
+ numerator += static_cast<RealType>( (*mapIt).second.m_TargetComplement );
+ denominator += static_cast<RealType>( (*mapIt).second.m_Target );
+ }
+ return ( numerator / denominator );
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetFalseNegativeError( LabelType label )
+{
+ MapIterator mapIt = this->m_LabelSetMeasures.find( label );
+ if( mapIt == this->m_LabelSetMeasures.end() )
+ {
+ itkWarningMacro( "Label " << label << " not found." );
+ return 0.0;
+ }
+ RealType value =
+ static_cast<RealType>( (*mapIt).second.m_TargetComplement ) /
+ static_cast<RealType>( (*mapIt).second.m_Target );
+ return value;
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetFalsePositiveError()
+{
+ RealType numerator = 0.0;
+ RealType denominator = 0.0;
+ for( MapIterator mapIt = this->m_LabelSetMeasures.begin();
+ mapIt != this->m_LabelSetMeasures.end(); ++mapIt )
+ {
+ // Do not include the background in the final value.
+ if( (*mapIt).first == NumericTraits<LabelType>::Zero )
+ {
+ continue;
+ }
+ numerator += static_cast<RealType>( (*mapIt).second.m_SourceComplement );
+ denominator += static_cast<RealType>( (*mapIt).second.m_Source );
+ }
+ return ( numerator / denominator );
+}
+
+template<class TLabelImage>
+typename LabelOverlapMeasuresImageFilter<TLabelImage>::RealType
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::GetFalsePositiveError( LabelType label )
+{
+ MapIterator mapIt = this->m_LabelSetMeasures.find( label );
+ if( mapIt == this->m_LabelSetMeasures.end() )
+ {
+ itkWarningMacro( "Label " << label << " not found." );
+ return 0.0;
+ }
+ RealType value =
+ static_cast<RealType>( (*mapIt).second.m_SourceComplement ) /
+ static_cast<RealType>( (*mapIt).second.m_Source );
+ return value;
+}
+
+template<class TLabelImage>
+void
+LabelOverlapMeasuresImageFilter<TLabelImage>
+::PrintSelf( std::ostream& os, Indent indent ) const
+{
+ Superclass::PrintSelf( os, indent );
+
+}
+
+
+}// end namespace itk
+#endif
motion_mask=$mask_dir/mm_$phase_nb.mhd
vf_result=$vf_dir/vf_${ref_phase_nb}_$phase_nb.mhd
clitkCombineImage -i $vf_in -j $vf_out -m $motion_mask -o $vf_result
+ abort_on_error registration $? clean_up_registration
+
clitkZeroVF -i $vf_in -o vf_zero.mhd
+ abort_on_error registration $? clean_up_registration
+
+ patient_mask=$mask_dir/patient_mask_$phase_nb.mhd
clitkCombineImage -i $vf_result -j vf_zero.mhd -m $patient_mask -o $vf_result
+ abort_on_error registration $? clean_up_registration
+
rm vf_zero.*
# save for later...
# create (zero) vf from reference to reference
clitkZeroVF -i $vf_ref -o $vf_dir/vf_${ref_phase_nb}_${ref_phase_nb}.mhd
+ abort_on_error registration $? clean_up_registration
# create 4D vf
create_mhd_4D_pattern.sh $vf_dir/vf_${ref_phase_nb}_
midp=$midp_dir/midp_$phase_nb.mhd
# average the vf's from reference phase to phase
clitkAverageTemporalDimension -i $vf_dir/vf_${ref_phase_nb}_4D.mhd -o $vf_midp
+ abort_on_error midp $? clean_up_midp
+
# invert the vf (why?)
clitkInvertVF -i $vf_midp -o $vf_midp
+ abort_on_error midp $? clean_up_midp
+
# create the midp by warping the reference phase with the reference vf
clitkWarpImage -i $ref_phase_file -o $midp --vf=$vf_midp -s 1
+ abort_on_error midp $? clean_up_midp
ref_vf_midp=$vf_midp
ref_midp=$midp
clitkImageConvert -i $ref_midp -o $ref_midp -t float
+ abort_on_error midp $? clean_up_midp
########### calculate the midp wrt the other phases
for i in $( seq 0 $((${#phase_files[@]} - 1))); do
# calculate vf from phase to midp, using the vf from reference phase to midp (-i)
# and the vf from reference phase to phase (-j)
clitkComposeVF -i $ref_vf_midp -j $vf_dir/vf_$ref_phase_nb\_$phase_nb.mhd -o $vf_midp
+ abort_on_error midp $? clean_up_midp
+
clitkWarpImage -i $phase_file -o $midp --vf=$vf_midp -s 1
+ abort_on_error midp $? clean_up_midp
+
clitkImageConvert -i $midp -o $midp -t float
+ abort_on_error midp $? clean_up_midp
fi
done
create_mhd_4D_pattern.sh $midp_dir/midp_
+
echo "Calculating midp_avg.mhd..."
clitkAverageTemporalDimension -i $midp_dir/midp_4D.mhd -o $midp_dir/midp_avg.mhd
+ abort_on_error midp $? clean_up_midp
+
echo "Calculating midp_med.mhd..."
clitkMedianTemporalDimension -i $midp_dir/midp_4D.mhd -o $midp_dir/midp_med.mhd
+ abort_on_error midp $? clean_up_midp
# clean-up
rm $midp_dir/vf_*
-#! /bin/bash -x
-
+#! /bin/bash
+
###############################################################################
#
# FILE: create_midP-2.0.sh
{
echo "$phase_file -> Extracting patient..."
clitkExtractPatient -i $phase_file -o $mask_dir_tmp/patient_mask_$phase_nb.mhd --noAutoCrop -a $afdb_file $ExtractPatientExtra
+# abort_on_error clitkExtractPatient $?
+
clitkSetBackground -i $phase_file -o $mask_dir_tmp/patient_$phase_nb.mhd --mask $mask_dir_tmp/patient_mask_$phase_nb.mhd --outsideValue -1000
+# abort_on_error clitkSetBackground $?
}
extract_bones()
extract_lungs()
{
- echo "$phase_file -> Extracting lungs..."
- clitkExtractLung -i $phase_file -o $mask_dir_tmp/lungs_$phase_nb.mhd -a $afdb_file --noAutoCrop
+ echo "$phase_file -> Extracting lungs..."
+ clitkExtractLung -i $phase_file -o $mask_dir_tmp/lungs_$phase_nb.mhd -a $afdb_file --noAutoCrop --doNotSeparateLungs
}
+
+
resample()
{
echo "$phase_file -> Resampling..."
clitkResampleImage -i $mask_dir_tmp/patient_$phase_nb.mhd -o $mask_dir_tmp/patient_$phase_nb.mhd --spacing $resample_spacing --interp $resample_algo
-
clitkResampleImage -i $mask_dir_tmp/lungs_$phase_nb.mhd -o $mask_dir_tmp/lungs_$phase_nb.mhd --like $mask_dir_tmp/patient_$phase_nb.mhd
}
mkdir -p $mask_log_dir
# multi-threaded pre-processing for motion mask calcs
+ pids=( )
for i in $( seq 0 $((${#phase_nbs[@]} - 1))); do
phase_nb=${phase_nbs[$i]}
phase_file=${phase_files[$i]}
check_threads $MAX_THREADS
mm_preprocessing &
+ pids=( "${pids[@]}" "$!" )
+ done
+
+ wait_pids ${pids[@]}
+ for ret in $ret_codes; do
+ abort_on_error mm_preprocessing $ret clean_up_masks
done
# single-threaded motion mask calc
check_threads 1
echo "$phase_file -> Computing motion mask..."
compute_motion_mask > $mask_log_dir/motion_mask_$phase_file.log
+ abort_on_error compute_motion_mask $? clean_up_masks
done
# multi-threaded post-processing of motion mask calcs
+ pids=( )
for i in $( seq 0 $((${#phase_nbs[@]} - 1))); do
phase_nb=${phase_nbs[$i]}
phase_file=${phase_files[$i]}
check_threads $MAX_THREADS
mm_postprocessing &
+ pids=( "${pids[@]}" "$!" )
done
+
+ wait_pids ${pids[@]}
+ for ret in $ret_codes; do
+ abort_on_error mm_postprocessing $ret clean_up_masks
+ done
+
# rename tmp mask directory after mask creation
check_threads 1
-#! /bin/sh +x
+#! /bin/sh -x
###############################################################################
#
#
###############################################################################
+#
+# check return value passed and abort if it represents an error (ie, ret != 0)
+# optionally, a function can be passed as a 3rd parameter, to be called just
+# before exiting. this is useful for cleaning up, for example.
+#
+abort_on_error()
+{
+ if [ $2 != 0 ]; then
+ echo Aborted at $1 with code $2
+ if [ $# = 3 ]; then
+ eval $3
+ fi
+ exit $2
+ fi
+}
+
+#
+# wait for all processes in the list and return their exit codes
+# in the ret_codes variable.
+#
+# OBS: this function must always be called in the same shell
+# that launched the processes.
+#
+wait_pids()
+{
+ local pids=$*
+ local ret=
+ local rets=
+# echo PIDS: $pids
+ for p in $pids; do
+# echo waiting $p
+ wait $p > /dev/null 2> /dev/null
+ ret=$?
+ if [ ret != 127 ]; then
+ rets=$rets" "$ret
+ else
+ rets=$rets" "0
+ fi
+
+ done
+
+ ret_codes=$rets
+}
+
+#
+# clean-up functions for maks, registration, and midp
+#
+clean_up_masks()
+{
+ rm -fr $mask_dir_tmp
+}
+
+clean_up_midp()
+{
+ rm -fr $midp_dir
+}
+
+clean_up_registration()
+{
+ rm -fr $vf_dir
+ rm -fr $output_dir
+}
+
+
+#
# block execution untill the number of threads (jobs) launched by the
# current process is below the given number of threads.
MAX_THREADS=2
#
###############################################################################
+source `dirname $0`/midp_common.sh
+
################# BLUTDIR #####################
registration_blutdir()
blutdir_params="--levels $nb_levels --metric $metric --optimizer $optimizer --samples $nb_samples --spacing $spacing,$spacing,$spacing --bins $hist_bins --maxIt $nb_iter --interp $interpolator --verbose"
cmd="clitkBLUTDIR -r $reference -t $target -m $mask_ref --targetMask $mask_targ --vf $vf -o $result $blutdir_params"
$cmd > $log
+
+ abort_on_error registration_blutdir $? clean_up_registration
}
################# ELASTIX #####################
# image registration
cmd="elastix -f $reference -m $target -fMask $mask_ref -mMask $mask_targ -out $result_dir -p params_BSpline_${suffix}.txt"
$cmd > /dev/null
+ abort_on_error registration_elastix $? clean_up_registration
# generate vector field
cmd="transformix -tp $result_dir/TransformParameters.0.txt -out $vf_dir -def all"
$cmd > /dev/null
+ abort_on_error registration_elastix $? clean_up_registration
# post-processing
mv $vf_dir/deformationField.mhd $vf
//--------------------------------------------------------------------
clitk::AnatomicalFeatureDatabase::AnatomicalFeatureDatabase()
{
- SetFilename("default.afdb");
+ SetFilename("default.afdb");
+ SetPath("./");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+clitk::AnatomicalFeatureDatabase::Pointer clitk::AnatomicalFeatureDatabase::New(std::string filename)
+{
+ Pointer a = AnatomicalFeatureDatabase::New();
+ a->SetFilename(filename);
+ a->Load();
+ return a;
}
//--------------------------------------------------------------------
//--------------------------------------------------------------------
+//--------------------------------------------------------------------
+std::string clitk::AnatomicalFeatureDatabase::GetTagValue(std::string tag)
+{
+ if (!TagExist(tag)) {
+ clitkExceptionMacro("Could not find the tag <" << tag << "> in the DB");
+ return "";
+ }
+ std::string s = m_MapOfTag[tag];
+ return s;
+}
+//--------------------------------------------------------------------
+
+
//--------------------------------------------------------------------
void clitk::AnatomicalFeatureDatabase::GetPoint3D(std::string tag, PointType3D & p)
{
class AnatomicalFeatureDatabase: public itk::Object
{
public:
- AnatomicalFeatureDatabase();
-
+ // typedef
typedef std::string TagType;
+ // New macro
+ typedef itk::Object Superclass;
+ typedef AnatomicalFeatureDatabase Self;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+ itkNewMacro(Self);
+ static Pointer New(std::string filename);
+
// Set/Get filename
itkSetMacro(Filename, std::string);
itkGetConstMacro(Filename, std::string);
// Read and write DB
void Write();
void Load();
+ itkGetConstMacro(Path, std::string);
+ itkSetMacro(Path, std::string);
// Set Get landmarks
typedef itk::Point<double,3> PointType3D;
void GetPoint3D(TagType tag, PointType3D & p);
double GetPoint3D(std::string tag, int dim);
bool TagExist(std::string tag);
-
+ std::string GetTagValue(std::string tag);
+
// Set Get image
void SetImageFilename(TagType tag, std::string f);
template<class ImageType>
void RemoveTag(TagType tag);
protected:
+ AnatomicalFeatureDatabase();
+ ~AnatomicalFeatureDatabase() {}
+
std::string m_Filename;
+ std::string m_Path;
typedef itk::ImageBase<3> ImageBaseType;
typedef std::map<TagType, std::string> MapTagType;
- typedef std::map<TagType, ImageBaseType*> MapTagImageType;
+ typedef std::map<TagType, ImageBaseType*> MapTagImageType;
MapTagType m_MapOfTag;
MapTagImageType m_MapOfImage;
}; // end class
//--------------------------------------------------------------------
+#define NewAFDB(NAME, FILE) \
+ clitk::AnatomicalFeatureDatabase::Pointer NAME = clitk::AnatomicalFeatureDatabase::New(FILE);
+
#include "clitkAnatomicalFeatureDatabase.txx"
} // end namespace clitk
else {
std::string s = m_MapOfTag[tag];
// Read the file
- image = readImage<ImageType>(s);
+ image = readImage<ImageType>(GetPath()+"/"+s);
// I add a reference count because the cache is not a smartpointer
image->SetReferenceCount(image->GetReferenceCount()+1);
// Insert into the cache
filter->SetArgsInfo(args_info);
- try {
- filter->Update();
- } catch(std::runtime_error e) {
- std::cout << e.what() << std::endl;
- }
+ CLITK_TRY_CATCH_EXIT(filter->Update());
return EXIT_SUCCESS;
} // This is the end, my friend
filter->SetArgsInfo(args_info);
- try {
- filter->Update();
- } catch(std::runtime_error e) {
- std::cerr << e.what() << std::endl;
- }
+ CLITK_TRY_CATCH_EXIT(filter->Update());
return EXIT_SUCCESS;
} // This is the end, my friend
option "doNotFillHoles" - "Do not fill holes if set" flag on
option "dir" d "Directions (axes) to perform filling (defaults to 2,1,0)" int multiple no
+section "Step 7 : lung separation (labelling)"
+option "doNotSeparateLungs" - "Do not separate lungs if set" flag off
+
option "noAutoCrop" - "If set : do no crop final mask to BoundingBox" flag off
itkGetConstMacro(FillHolesFlag, bool);
itkBooleanMacro(FillHolesFlag);
+ // Separate lungs
+ itkSetMacro(SeparateLungsFlag, bool);
+ itkGetConstMacro(SeparateLungsFlag, bool);
+ itkBooleanMacro(SeparateLungsFlag);
+
// Step Auto Crop
itkSetMacro(AutoCrop, bool);
itkGetConstMacro(AutoCrop, bool);
bool m_FillHolesFlag;
InputImageSizeType m_FillHolesDirections;
+ bool m_SeparateLungsFlag;
+
// Main functions
virtual void GenerateOutputInformation();
virtual void GenerateInputRequestedRegion();
StopCurrentStep<MaskImageType>(working_mask);
}
- //--------------------------------------------------------------------
- //--------------------------------------------------------------------
- StartNewStep("Separate Left/Right lungs");
- PrintMemory(GetVerboseMemoryFlag(), "Before Separate");
- // Initial label
- working_mask = Labelize<MaskImageType>(working_mask,
- GetBackgroundValue(),
- false,
- GetMinimalComponentSize());
-
- PrintMemory(GetVerboseMemoryFlag(), "After Labelize");
-
- // Count the labels
- typedef itk::StatisticsImageFilter<MaskImageType> StatisticsImageFilterType;
- typename StatisticsImageFilterType::Pointer statisticsImageFilter=StatisticsImageFilterType::New();
- statisticsImageFilter->SetInput(working_mask);
- statisticsImageFilter->Update();
- unsigned int initialNumberOfLabels = statisticsImageFilter->GetMaximum();
- working_mask = statisticsImageFilter->GetOutput();
+ if (GetSeparateLungsFlag()) {
+ //--------------------------------------------------------------------
+ //--------------------------------------------------------------------
+ StartNewStep("Separate Left/Right lungs");
+ PrintMemory(GetVerboseMemoryFlag(), "Before Separate");
+ // Initial label
+ working_mask = Labelize<MaskImageType>(working_mask,
+ GetBackgroundValue(),
+ false,
+ GetMinimalComponentSize());
+
+ PrintMemory(GetVerboseMemoryFlag(), "After Labelize");
+
+ // Count the labels
+ typedef itk::StatisticsImageFilter<MaskImageType> StatisticsImageFilterType;
+ typename StatisticsImageFilterType::Pointer statisticsImageFilter=StatisticsImageFilterType::New();
+ statisticsImageFilter->SetInput(working_mask);
+ statisticsImageFilter->Update();
+ unsigned int initialNumberOfLabels = statisticsImageFilter->GetMaximum();
+ working_mask = statisticsImageFilter->GetOutput();
+
+ PrintMemory(GetVerboseMemoryFlag(), "After count label");
- PrintMemory(GetVerboseMemoryFlag(), "After count label");
-
- // Decompose the first label
- if (initialNumberOfLabels<2) {
- // Structuring element radius
- typename ImageType::SizeType radius;
- for (unsigned int i=0;i<Dim;i++) radius[i]=1;
- typedef clitk::DecomposeAndReconstructImageFilter<MaskImageType,MaskImageType> DecomposeAndReconstructFilterType;
- typename DecomposeAndReconstructFilterType::Pointer decomposeAndReconstructFilter=DecomposeAndReconstructFilterType::New();
- decomposeAndReconstructFilter->SetInput(working_mask);
- decomposeAndReconstructFilter->SetVerbose(false);
- decomposeAndReconstructFilter->SetRadius(radius);
- decomposeAndReconstructFilter->SetMaximumNumberOfLabels(2);
- decomposeAndReconstructFilter->SetMinimumObjectSize(this->GetMinimalComponentSize());
- decomposeAndReconstructFilter->SetMinimumNumberOfIterations(1);
- decomposeAndReconstructFilter->SetBackgroundValue(this->GetBackgroundValue());
- decomposeAndReconstructFilter->SetForegroundValue(this->GetForegroundValue());
- decomposeAndReconstructFilter->SetFullyConnected(true);
- decomposeAndReconstructFilter->SetNumberOfNewLabels(1);
- decomposeAndReconstructFilter->Update();
- working_mask = decomposeAndReconstructFilter->GetOutput();
+ // Decompose the first label
+ if (initialNumberOfLabels<2) {
+ // Structuring element radius
+ typename ImageType::SizeType radius;
+ for (unsigned int i=0;i<Dim;i++) radius[i]=1;
+ typedef clitk::DecomposeAndReconstructImageFilter<MaskImageType,MaskImageType> DecomposeAndReconstructFilterType;
+ typename DecomposeAndReconstructFilterType::Pointer decomposeAndReconstructFilter=DecomposeAndReconstructFilterType::New();
+ decomposeAndReconstructFilter->SetInput(working_mask);
+ decomposeAndReconstructFilter->SetVerbose(false);
+ decomposeAndReconstructFilter->SetRadius(radius);
+ decomposeAndReconstructFilter->SetMaximumNumberOfLabels(2);
+ decomposeAndReconstructFilter->SetMinimumObjectSize(this->GetMinimalComponentSize());
+ decomposeAndReconstructFilter->SetMinimumNumberOfIterations(1);
+ decomposeAndReconstructFilter->SetBackgroundValue(this->GetBackgroundValue());
+ decomposeAndReconstructFilter->SetForegroundValue(this->GetForegroundValue());
+ decomposeAndReconstructFilter->SetFullyConnected(true);
+ decomposeAndReconstructFilter->SetNumberOfNewLabels(1);
+ decomposeAndReconstructFilter->Update();
+ working_mask = decomposeAndReconstructFilter->GetOutput();
+ }
+ PrintMemory(GetVerboseMemoryFlag(), "After decomposeAndReconstructFilter");
}
- PrintMemory(GetVerboseMemoryFlag(), "After decomposeAndReconstructFilter");
// Retain labels ('1' is largset lung, so right. '2' is left)
typedef itk::ThresholdImageFilter<MaskImageType> ThresholdImageFilterType;
f->SetFillHolesFlag(false);
else
f->SetFillHolesFlag(true);
+
+ if (mArgsInfo.doNotSeparateLungs_given)
+ f->SetSeparateLungsFlag(false);
+ else
+ f->SetSeparateLungsFlag(true);
}
//--------------------------------------------------------------------
Aorta = clitk::Dilate<MaskImageType>(Aorta, radius, GetBackgroundValue(), GetForegroundValue(), false);
// Now, insert this image in the AFDB (but do not store on disk)
- GetAFDB()->template SetImage<MaskImageType>("Aorta_Dilated_Anteriorly", "bidon",
- Aorta, false);
- /*
- // Not Post to
- m_Working_Support =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_Working_Support, Aorta, 2,
- GetFuzzyThreshold("3A", "Aorta"),
- "NotPostTo", true,
- Aorta->GetSpacing()[0], false, false);
-
- */
+ GetAFDB()->template SetImage<MaskImageType>("Aorta_Dilated_Anteriorly", "seg/Aorta_Dilated_Anteriorly.mha", Aorta, false);
+ writeImage<MaskImageType>(Aorta, "seg/Aorta_Dilated_Anteriorly.mha");
+ GetAFDB()->Write();
-
StopCurrentStep<MaskImageType>(m_Working_Support);
m_ListOfStations["3A"] = m_Working_Support;
}
template <class TImageType>
void
clitk::ExtractLymphStationsFilter<TImageType>::
-ExtractStation_4RL() {
- if ((!CheckForStation("4R")) && (!CheckForStation("4L"))) return;
-
- StartNewStep("Stations 4RL");
+ExtractStation_4R() {
+ if (!CheckForStation("4R")) return;
+ StartNewStep("Stations 4R");
StartSubStep();
// Get the current support
- StartNewStep("[Station 4RL] Get the current 4RL suppport");
+ StartNewStep("[Station 4R] Get the current 4RL suppport");
m_ListOfStations["4R"] = m_ListOfSupports["S4R"];
- m_ListOfStations["4L"] = m_ListOfSupports["S4L"];
StopCurrentStep<MaskImageType>(m_ListOfStations["4R"]);
// Generic RelativePosition processes
m_ListOfStations["4R"] = this->ApplyRelativePositionList("Station_4R", m_ListOfStations["4R"]);
+
+ // Store image filenames into AFDB
+ WriteImageStation("4R");
+ StopSubStep();
+ ComputeOverlapWithRef("4R");
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class TImageType>
+void
+clitk::ExtractLymphStationsFilter<TImageType>::
+ExtractStation_4L() {
+ if (!CheckForStation("4L")) return;
+ StartNewStep("Stations 4L");
+ StartSubStep();
+
+ // Get the current support
+ StartNewStep("[Station 4L] Get the current 4RL suppport");
+ m_ListOfStations["4L"] = m_ListOfSupports["S4L"];
+ StopCurrentStep<MaskImageType>(m_ListOfStations["4L"]);
+
+ // Generic RelativePosition processes
m_ListOfStations["4L"] = this->ApplyRelativePositionList("Station_4L", m_ListOfStations["4L"]);
// Separation Ant/Post
StopCurrentStep<MaskImageType>(m_ListOfStations["4L"]);
// Store image filenames into AFDB
- writeImage<MaskImageType>(m_ListOfStations["4R"], "seg/Station4R.mhd");
- writeImage<MaskImageType>(m_ListOfStations["4L"], "seg/Station4L.mhd");
- GetAFDB()->SetImageFilename("Station4R", "seg/Station4R.mhd");
- GetAFDB()->SetImageFilename("Station4L", "seg/Station4L.mhd");
- WriteAFDB();
+ WriteImageStation("4L");
StopSubStep();
-
+ ComputeOverlapWithRef("4L");
}
//--------------------------------------------------------------------
// Get initial Mediastinum
m_Working_Support = m_Mediastinum = this->GetAFDB()->template GetImage<MaskImageType>("Mediastinum", true);
- // Consider sup/inf to Carina
- double m_CarinaZ = FindCarina();
- MaskImagePointer m_Support_Superior_to_Carina =
- clitk::CropImageRemoveLowerThan<MaskImageType>(m_Working_Support, 2,
- m_CarinaZ, true, GetBackgroundValue());
- MaskImagePointer m_Support_Inferior_to_Carina =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(m_Working_Support, 2,
- m_CarinaZ, true, GetBackgroundValue());
- m_ListOfSupports["Support_Superior_to_Carina"] = m_Support_Superior_to_Carina;
- m_ListOfSupports["Support_Inferior_to_Carina"] = m_Support_Inferior_to_Carina;
- writeImage<MaskImageType>(m_Support_Inferior_to_Carina, "seg/Support_Inf_Carina.mha");
- this->GetAFDB()->SetImageFilename("Support_Inf_Carina", "seg/Support_Inf_Carina.mha");
- writeImage<MaskImageType>(m_Support_Superior_to_Carina, "seg/Support_Sup_Carina.mha");
- this->GetAFDB()->SetImageFilename("Support_Sup_Carina", "seg/Support_Sup_Carina.mha");
+ // Remove some computed structures
+ this->GetAFDB()->RemoveTag("CarinaZ");
+ this->GetAFDB()->RemoveTag("ApexOfTheChestZ");
+
+ // Superior and inferior limits.
+ Support_SI_Limit("inferior", "Sup_to_Carina", "inferior", "Carina", 0);
+ Support_SI_Limit("superior", "Inf_to_Carina", "inferior", "Carina", m_Working_Support->GetSpacing()[2]);
+
+ // Initialise all others supports
+ // m_ListOfSupports["S1R"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S1L"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S2R"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S2L"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S3A"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S3P"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S4R"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S4L"] = m_ListOfSupports["Sup_to_Carina"];
+ // m_ListOfSupports["S5"] = m_Mediastinum; // Not above Carina
+ // m_ListOfSupports["S6"] = m_Mediastinum; // Not above Carina
+
+ // Read all support limits in a file and apply them
+ ReadSupportLimits(GetSupportLimitsFilename());
+ for(int i=0; i<m_ListOfSupportLimits.size(); i++) {
+ SupportLimitsType s = m_ListOfSupportLimits[i];
+ Support_SI_Limit(s.station_limit, s.station, s.structure_limit,
+ s.structure, s.offset*m_Working_Support->GetSpacing()[2]);
+ }
// S1RL
- Support_SupInf_S1RL();
Support_LeftRight_S1R_S1L();
// S2RL
- Support_SupInf_S2R_S2L();
Support_LeftRight_S2R_S2L();
// S4RL
- Support_SupInf_S4R_S4L();
Support_LeftRight_S4R_S4L();
// Post limits of S1,S2,S4
Support_Post_S1S2S4();
- // S3AP
- Support_S3P();
- Support_S3A();
-
- // S5, S6
- Support_S5();
- Support_S6();
+ // S3P
+ StartNewStep("[Support] Ant limits of S3P with trachea");
+ m_ListOfSupports["S3P"] = LimitsWithTrachea(m_ListOfSupports["S3P"], 1, 0, 10);
+
+ // S3A
+ StartNewStep("[Support] Ant limits of S3A with trachea");
+ m_ListOfSupports["S3A"] = LimitsWithTrachea(m_ListOfSupports["S3A"], 1, 0, -10);
+ // I will do it later
// Below Carina S7,8,9,10
- m_ListOfSupports["S7"] = clitk::Clone<MaskImageType>(m_Support_Inferior_to_Carina);
- m_ListOfSupports["S8"] = clitk::Clone<MaskImageType>(m_Support_Inferior_to_Carina);
- m_ListOfSupports["S9"] = clitk::Clone<MaskImageType>(m_Support_Inferior_to_Carina);
- m_ListOfSupports["S10"] = clitk::Clone<MaskImageType>(m_Support_Inferior_to_Carina);
- m_ListOfSupports["S11"] = clitk::Clone<MaskImageType>(m_Support_Inferior_to_Carina);
+ m_ListOfSupports["S7"] = clitk::Clone<MaskImageType>(m_ListOfSupports["Inf_to_Carina"]);
+ m_ListOfSupports["S8"] = clitk::Clone<MaskImageType>(m_ListOfSupports["Inf_to_Carina"]);
+ m_ListOfSupports["S9"] = clitk::Clone<MaskImageType>(m_ListOfSupports["Inf_to_Carina"]);
+ m_ListOfSupports["S10"] = clitk::Clone<MaskImageType>(m_ListOfSupports["Inf_to_Carina"]);
+ m_ListOfSupports["S11"] = clitk::Clone<MaskImageType>(m_ListOfSupports["Inf_to_Carina"]);
// Store image filenames into AFDB
- writeImage<MaskImageType>(m_ListOfSupports["S1R"], "seg/Support_S1R.mha");
- this->GetAFDB()->SetImageFilename("Support_S1R", "seg/Support_S1R.mha");
- writeImage<MaskImageType>(m_ListOfSupports["S1L"], "seg/Support_S1L.mha");
- this->GetAFDB()->SetImageFilename("Support_S1L", "seg/Support_S1L.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S2L"], "seg/Support_S2L.mha");
- this->GetAFDB()->SetImageFilename("Support_S2L", "seg/Support_S2L.mha");
- writeImage<MaskImageType>(m_ListOfSupports["S2R"], "seg/Support_S2R.mha");
- this->GetAFDB()->SetImageFilename("Support_S2R", "seg/Support_S2R.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S3P"], "seg/Support_S3P.mha");
- this->GetAFDB()->SetImageFilename("Support_S3P", "seg/Support_S3P.mha");
- writeImage<MaskImageType>(m_ListOfSupports["S3A"], "seg/Support_S3A.mha");
- this->GetAFDB()->SetImageFilename("Support_S3A", "seg/Support_S3A.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S4L"], "seg/Support_S4L.mha");
- this->GetAFDB()->SetImageFilename("Support_S4L", "seg/Support_S4L.mha");
- writeImage<MaskImageType>(m_ListOfSupports["S4R"], "seg/Support_S4R.mha");
- this->GetAFDB()->SetImageFilename("Support_S4R", "seg/Support_S4R.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S5"], "seg/Support_S5.mha");
- this->GetAFDB()->SetImageFilename("Support_S5", "seg/Support_S5.mha");
- writeImage<MaskImageType>(m_ListOfSupports["S6"], "seg/Support_S6.mha");
- this->GetAFDB()->SetImageFilename("Support_S6", "seg/Support_S6.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S7"], "seg/Support_S7.mha");
- this->GetAFDB()->SetImageFilename("Support_S7", "seg/Support_S7.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S8"], "seg/Support_S8.mha");
- this->GetAFDB()->SetImageFilename("Support_S8", "seg/Support_S8.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S9"], "seg/Support_S9.mha");
- this->GetAFDB()->SetImageFilename("Support_S9", "seg/Support_S9.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S10"], "seg/Support_S10.mha");
- this->GetAFDB()->SetImageFilename("Support_S10", "seg/Support_S10.mha");
-
- writeImage<MaskImageType>(m_ListOfSupports["S11"], "seg/Support_S11.mha");
- this->GetAFDB()->SetImageFilename("Support_S11", "seg/Support_S11.mha");
+ WriteImageSupport("S1R"); WriteImageSupport("S1L");
+ WriteImageSupport("S2R"); WriteImageSupport("S2L");
+ WriteImageSupport("S3A"); WriteImageSupport("S3P");
+ WriteImageSupport("S4R"); WriteImageSupport("S4L");
+ WriteImageSupport("S5");
+ WriteImageSupport("S6");
+ WriteImageSupport("S7");
+ WriteImageSupport("S8");
+ WriteImageSupport("S9");
+ WriteImageSupport("S10");
+ WriteImageSupport("S11");
WriteAFDB();
}
//--------------------------------------------------------------------
template <class ImageType>
void
clitk::ExtractLymphStationsFilter<ImageType>::
-Support_SupInf_S1RL()
+Support_SI_Limit(const std::string station_limit, const std::string station,
+ const std::string structure_limit, const std::string structure,
+ const double offset)
{
- // Step : S1RL
- StartNewStep("[Support] Sup-Inf S1RL");
- /*
- 2R: Upper border: apex of the right lung and pleural space, and in
- the midline, the upper border of the manubrium
-
- 2L: Upper border: apex of the left lung and pleural space, and in the
- midline, the upper border of the manubrium
+ if (!GetCheckSupportFlag())
+ StartNewStep("[Support] "+station_limit+" limit of "+station+" is "+structure_limit+" limit of "+structure);
- => apex / manubrium = up Sternum
- */
- m_Working_Support = m_ListOfSupports["Support_Superior_to_Carina"];
- MaskImagePointer Sternum = this->GetAFDB()->template GetImage <MaskImageType>("Sternum");
- MaskImagePointType p;
- p[0] = p[1] = p[2] = 0.0; // to avoid warning
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Sternum, GetBackgroundValue(), 2, false, p);
- MaskImagePointer S1RL =
- clitk::CropImageRemoveLowerThan<MaskImageType>(m_Working_Support, 2,
- p[2], true, GetBackgroundValue());
- m_Working_Support =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(m_Working_Support, 2,
- p[2], true, GetBackgroundValue());
- m_ListOfSupports["S1RL"] = S1RL;
+ // Check
+ if ((station_limit != "superior") && (station_limit != "inferior")) {
+ clitkExceptionMacro("Error station_limit must be 'inferior' or 'superior', not '"<< station_limit);
+ }
+ if ((structure_limit != "superior") && (structure_limit != "inferior")) {
+ clitkExceptionMacro("Error structure_limit must be 'inferior' or 'superior', not '"<< structure_limit);
+ }
+
+ // Get current support
+ if (m_ListOfSupports.find(station) == m_ListOfSupports.end()) {
+ // std::cerr << "Warning: support " << station << " not initialized" << std::endl;
+ m_ListOfSupports[station] = m_Mediastinum;
+ }
+ m_Working_Support = m_ListOfSupports[station];
+
+ // Get structure or structureZ
+ double z;
+ int found=0;
+ std::string file;
+
+ // Try to load structure and compute extrema point
+ if (this->GetAFDB()->TagExist(structure)) {
+ MaskImagePointer Structure = this->GetAFDB()->template GetImage <MaskImageType>(structure);
+ file = this->GetAFDB()->GetTagValue(structure);
+ MaskImagePointType p;
+ p[0] = p[1] = p[2] = 0.0; // to avoid warning
+ if (structure_limit == "superior")
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Structure, GetBackgroundValue(), 2, false, p);
+ else
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Structure, GetBackgroundValue(), 2, true, p);
+ z = p[2];
+ found=1;
+ }
+
+ // Try to load structureZ
+ if ((found==0) && (this->GetAFDB()->TagExist(structure+"Z"))) {
+ z = this->GetAFDB()->GetDouble(structure+"Z");
+ found=2;
+ }
+
+ // Try to load structurePoint
+ if ((found==0) && (this->GetAFDB()->TagExist(structure+"Point"))) {
+ MaskImagePointType p;
+ this->GetAFDB()->GetPoint3D(structure+"Point", p);
+ z = p[2];
+ found=3;
+ }
+
+ // Try to see if it is an already computed support
+ if (found==0) {
+ if (m_ListOfSupports.find(structure) != m_ListOfSupports.end()) {
+ MaskImagePointer Structure = m_ListOfSupports[structure];
+ MaskImagePointType p;
+ if (structure_limit == "superior")
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Structure, GetBackgroundValue(), 2, false, p);
+ else
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Structure, GetBackgroundValue(), 2, true, p);
+ z = p[2];
+ found=4;
+ }
+ }
+
+ // Try special case : "FindApexOfTheChest"
+ if (structure == "FindApexOfTheChest") {
+ z = FindApexOfTheChest();
+ found=5;
+ }
+ if (structure == "FindInferiorBorderOfAorticArch") {
+ z = FindInferiorBorderOfAorticArch();
+ found=6;
+ }
+ if (structure == "FindSuperiorBorderOfAorticArch") {
+ z = FindSuperiorBorderOfAorticArch();
+ found=6;
+ }
+
+ // If we find anything
+ if (found == 0) {
+ std::cerr << "ERROR : I could not find " << structure << " nor " << structure << "Z nor "
+ << structure << "Point" << std::endl;
+ exit(EXIT_FAILURE);
+ }
+
+ // Apply offset
+ z += offset;
+
+ // Remove Lower or greater
+ if (station_limit == "inferior") {
+ m_Working_Support =
+ clitk::CropImageRemoveLowerThan<MaskImageType>(m_Working_Support, 2, z, true, GetBackgroundValue());
+ }
+ else {
+ m_Working_Support =
+ clitk::CropImageRemoveGreaterThan<MaskImageType>(m_Working_Support, 2, z, true, GetBackgroundValue());
+ }
+
+ // Check: if reference station is given, display information
+ if (GetCheckSupportFlag()) {
+ try {
+ MaskImagePointer Ref = this->GetAFDB()->template GetImage <MaskImageType>(station+"_Ref");
+ MaskImagePointType p_support;
+ MaskImagePointType p_ref;
+ if (station_limit == "superior") {
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Ref, GetBackgroundValue(), 2, false, p_ref);
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(m_Working_Support, GetBackgroundValue(), 2, false, p_support);
+ }
+ else {
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Ref, GetBackgroundValue(), 2, true, p_ref);
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(m_Working_Support, GetBackgroundValue(), 2, true, p_support);
+ }
+ std::ostringstream os;
+ os << "[Support] \t" << station << "\t" << station_limit << " "
+ << "Z = " << z << std::setprecision(2) << std::fixed
+ << "\tSupport = " << p_support[2]
+ << "\tRef = " << p_ref[2]
+ << "\tdiff = " << p_support[2]-p_ref[2] << "\t"
+ << structure << " " << structure_limit;
+ if (found==1) os << " (S "+file+")";
+ if (found==2) os << " (Z)";
+ if (found==3) os << " (P)";
+ if (found==4) os << " (p)";
+ if (found==5) os << " (Apex)";
+ if (found==6) os << " (AorticArch)";
+ StartNewStep(os.str());
+ } catch(clitk::ExceptionObject e) { }
+ }
+
+ // Set support
+ m_ListOfSupports[station] = m_Working_Support;
+ StopCurrentStep<MaskImageType>(m_Working_Support);
}
//--------------------------------------------------------------------
std::vector<ImagePointType> B;
// Search for centroid positions of trachea
MaskImagePointer Trachea = this->GetAFDB()->template GetImage <MaskImageType>("Trachea");
- MaskImagePointer S1RL = m_ListOfSupports["S1RL"];
+ MaskImagePointer S1RL = m_ListOfSupports["S1R"];
Trachea = clitk::ResizeImageLike<MaskImageType>(Trachea, S1RL, GetBackgroundValue());
std::vector<MaskSlicePointer> slices;
clitk::ExtractSlices<MaskImageType>(Trachea, 2, slices);
// Right part
clitk::SliceBySliceSetBackgroundFromLineSeparation<MaskImageType>(S1R, A, B,
- GetBackgroundValue(), 0, -10);
+ GetBackgroundValue(), 0, 10);
S1R = clitk::AutoCrop<MaskImageType>(S1R, GetBackgroundValue());
m_ListOfSupports["S1R"] = S1R;
// Left part
clitk::SliceBySliceSetBackgroundFromLineSeparation<MaskImageType>(S1L, A, B,
- GetBackgroundValue(), 0, 10);
+ GetBackgroundValue(), 0, -10);
S1L = clitk::AutoCrop<MaskImageType>(S1L, GetBackgroundValue());
m_ListOfSupports["S1L"] = S1L;
+ StopCurrentStep<MaskImageType>(m_ListOfSupports["S1L"]);
}
//--------------------------------------------------------------------
-//--------------------------------------------------------------------
-template <class ImageType>
-void
-clitk::ExtractLymphStationsFilter<ImageType>::
-Support_SupInf_S2R_S2L()
-{
- // Step : S2RL Sup-Inf limits
- /*
- 2R Lower border: intersection of caudal margin of innominate vein with
- the trachea
- 2L Lower border: superior border of the aortic arch
- */
- StartNewStep("[Support] Sup-Inf S2RL");
- m_Working_Support = m_ListOfSupports["Support_Superior_to_Carina"];
-
- // S2R Caudal Margin Of Left BrachiocephalicVein
- MaskImagePointer BrachioCephalicVein = this->GetAFDB()->template GetImage<MaskImageType>("BrachioCephalicVein");
- MaskImagePointType p;
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(BrachioCephalicVein, GetBackgroundValue(), 2, true, p);
-
- // I add slightly more than a slice --> NO !!
- double CaudalMarginOfLeftBrachiocephalicVeinZ=p[2];//+ 1.1*m_Working_Support->GetSpacing()[2];
-
- this->GetAFDB()->SetDouble("CaudalMarginOfLeftBrachiocephalicVeinZ", CaudalMarginOfLeftBrachiocephalicVeinZ);
- MaskImagePointer S2R =
- clitk::CropImageRemoveLowerThan<MaskImageType>(m_Working_Support, 2,
- CaudalMarginOfLeftBrachiocephalicVeinZ, true,
- GetBackgroundValue());
- // S2L : Top Of Aortic Arch
- MaskImagePointer Aorta = this->GetAFDB()->template GetImage<MaskImageType>("Aorta");
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Aorta, GetBackgroundValue(), 2, false, p);
-
- // Save the TopOfAorticArchZ
- this->GetAFDB()->SetDouble("TopOfAorticArchZ", p[2]);
-
- // I substract slightly more than a slice to respect delineation
- double TopOfAorticArchZ=p[2]- 1.1*m_Working_Support->GetSpacing()[2];
- // this->GetAFDB()->SetDouble("TopOfAorticArchZ", TopOfAorticArchZ);
-
- MaskImagePointer S2L =
- clitk::CropImageRemoveLowerThan<MaskImageType>(m_Working_Support, 2,
- TopOfAorticArchZ, true,
- GetBackgroundValue());
-
- /*
- // S2RL: Superior support, I use inferior part of S1RL
- MaskImagePointer S1L = m_ListOfSupports["S1L"];
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(S1L, GetBackgroundValue(), 2, true, p);
- DD(p);
- S2L =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(S2L, 2,
- p[2], true,
- GetBackgroundValue());
-
- MaskImagePointer S1R = m_ListOfSupports["S1R"];
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(S1R, GetBackgroundValue(), 2, true, p);
- DD(p);
- S2R =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(S2R, 2,
- p[2], true,
- GetBackgroundValue());
- */
-
- // Superior limits, use Sternum (but not strictly inf to S1RL
- MaskImagePointer Sternum = this->GetAFDB()->template GetImage <MaskImageType>("Sternum");
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Sternum, GetBackgroundValue(), 2, false, p);
- // Add one slice
- p[2] = p[2] + m_Working_Support->GetSpacing()[2];
- S2L =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(S2L, 2,
- p[2], true, GetBackgroundValue());
- S2R =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(S2R, 2,
- p[2], true, GetBackgroundValue());
-
- // The is the end
- m_ListOfSupports["S2L"] = S2L;
- m_ListOfSupports["S2R"] = S2R;
-}
-//--------------------------------------------------------------------
-
-
-
//--------------------------------------------------------------------
template <class ImageType>
void
MaskImagePointer S2L = m_ListOfSupports["S2L"];
S2R = LimitsWithTrachea(S2R, 0, 1, -10);
S2L = LimitsWithTrachea(S2L, 0, 1, 10);
+ S2R = clitk::AutoCrop<MaskImageType>(S2R, GetBackgroundValue());
+ S2L = clitk::AutoCrop<MaskImageType>(S2L, GetBackgroundValue());
m_ListOfSupports["S2R"] = S2R;
m_ListOfSupports["S2L"] = S2L;
this->GetAFDB()->template ReleaseImage<MaskImageType>("Trachea");
//--------------------------------------------------------------------
-//--------------------------------------------------------------------
-template <class ImageType>
-void
-clitk::ExtractLymphStationsFilter<ImageType>::
-Support_SupInf_S4R_S4L()
-{
- // ---------------------------------------------------------------------------
- /* Step : S4RL Sup-Inf
- - start at the end of 2R and 2L
- - stop ?
- - 4R
- Rod says : "The inferior border is at the lower border of the azygous vein."
- Rod says : difficulties
- (was : "ends at the upper lobe bronchus or where the right pulmonary artery
- crosses the midline of the mediastinum ")
- - 4L
- Rod says : "The lower border is to upper margin of the left main pulmonary artery."
- (was LLL bronchus)
- */
- StartNewStep("[Support] Sup-Inf limits of 4R/4L");
-
- // Start from the support
- MaskImagePointer S4RL = clitk::Clone<MaskImageType>(m_Working_Support);
- MaskImagePointer S4R = clitk::Clone<MaskImageType>(S4RL);
- MaskImagePointer S4L = clitk::Clone<MaskImageType>(S4RL);
-
- // Keep only what is lower than S2
- MaskImagePointer S2R = m_ListOfSupports["S2R"];
- MaskImagePointer S2L = m_ListOfSupports["S2L"];
- MaskImagePointType p;
- // Right part
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(S2R, GetBackgroundValue(),
- 2, true, p);
- S4R = clitk::CropImageRemoveGreaterThan<MaskImageType>(S4R, 2,
- p[2], true, GetBackgroundValue());
- // Left part
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(S2L, GetBackgroundValue(),
- 2, true, p);
- S4L = clitk::CropImageRemoveGreaterThan<MaskImageType>(S4L, 2,
- p[2], true, GetBackgroundValue());
-
- // Get AzygousVein and limit according to LowerBorderAzygousVein
- MaskImagePointer LowerBorderAzygousVein
- = this->GetAFDB()->template GetImage<MaskImageType>("LowerBorderAzygousVein");
- std::vector<MaskImagePointType> c;
- clitk::ComputeCentroids<MaskImageType>(LowerBorderAzygousVein, GetBackgroundValue(), c);
- S4R = clitk::CropImageRemoveLowerThan<MaskImageType>(S4R, 2,
- c[1][2], true, GetBackgroundValue());
- S4R = clitk::AutoCrop<MaskImageType>(S4R, GetBackgroundValue());
- m_ListOfSupports["S4R"] = S4R;
-
-
- // Limit according to LeftPulmonaryArtery
- MaskImagePointer LeftPulmonaryArtery
- = this->GetAFDB()->template GetImage<MaskImageType>("LeftPulmonaryArtery");
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(LeftPulmonaryArtery, GetBackgroundValue(),
- 2, false, p);
- S4L = clitk::CropImageRemoveLowerThan<MaskImageType>(S4L, 2,
- p[2], true, GetBackgroundValue());
- S4L = clitk::AutoCrop<MaskImageType>(S4L, GetBackgroundValue());
- m_ListOfSupports["S4L"] = S4L;
-}
-//--------------------------------------------------------------------
-
-
//--------------------------------------------------------------------
template <class ImageType>
void
MaskImagePointer S2L = m_ListOfSupports["S2L"];
MaskImagePointer S4R = m_ListOfSupports["S4R"];
MaskImagePointer S4L = m_ListOfSupports["S4L"];
- S1L = LimitsWithTrachea(S1L, 1, 0, -10, m_ApexOfTheChest);
- S1R = LimitsWithTrachea(S1R, 1, 0, -10, m_ApexOfTheChest);
- S2R = LimitsWithTrachea(S2R, 1, 0, -10, m_ApexOfTheChest);
- S2L = LimitsWithTrachea(S2L, 1, 0, -10, m_ApexOfTheChest);
- S4R = LimitsWithTrachea(S4R, 1, 0, -10, m_ApexOfTheChest);
- S4L = LimitsWithTrachea(S4L, 1, 0, -10, m_ApexOfTheChest);
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class ImageType>
-void
-clitk::ExtractLymphStationsFilter<ImageType>::
-Support_S3P()
-{
- StartNewStep("[Support] Ant limits of S3P and Post limits of S1RL, S2RL, S4RL");
-
- // Initial S3P support
- MaskImagePointer S3P = clitk::Clone<MaskImageType>(m_ListOfSupports["Support_Superior_to_Carina"]);
-
- // Stop at Lung Apex
- double m_ApexOfTheChest = FindApexOfTheChest();
- S3P =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(S3P, 2,
- m_ApexOfTheChest, true,
- GetBackgroundValue());
- // Ant limits with Trachea
- S3P = LimitsWithTrachea(S3P, 1, 0, 10);
- m_ListOfSupports["S3P"] = S3P;
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class ImageType>
-void
-clitk::ExtractLymphStationsFilter<ImageType>::
-Support_S3A()
-{
- StartNewStep("[Support] Sup-Inf and Post limits for S3A");
-
- // Initial S3A support
- MaskImagePointer S3A = clitk::Clone<MaskImageType>(m_ListOfSupports["Support_Superior_to_Carina"]);
-
- // Stop at Lung Apex or like S2/S1 (upper border Sternum - manubrium) ?
-
- //double m_ApexOfTheChest = FindApexOfTheChest();
-
- MaskImagePointer Sternum = this->GetAFDB()->template GetImage <MaskImageType>("Sternum");
- MaskImagePointType p;
- p[0] = p[1] = p[2] = 0.0; // to avoid warning
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(Sternum, GetBackgroundValue(), 2, false, p);
- p[2] += Sternum->GetSpacing()[2]; // we add one slice to stop 3A at the same slice than Sternum stop
- S3A =
- clitk::CropImageRemoveGreaterThan<MaskImageType>(S3A, 2,
- //m_ApexOfTheChest
- p[2], true,
- GetBackgroundValue());
- // Ant limits with Trachea
- S3A = LimitsWithTrachea(S3A, 1, 0, -10);
- m_ListOfSupports["S3A"] = S3A;
+ m_ListOfSupports["S1R"] = LimitsWithTrachea(S1L, 1, 0, -10, m_ApexOfTheChest);
+ m_ListOfSupports["S1L"] = LimitsWithTrachea(S1R, 1, 0, -10, m_ApexOfTheChest);
+ m_ListOfSupports["S2R"] = LimitsWithTrachea(S2R, 1, 0, -10, m_ApexOfTheChest);
+ m_ListOfSupports["S2L"] = LimitsWithTrachea(S2L, 1, 0, -10, m_ApexOfTheChest);
+ m_ListOfSupports["S4R"] = LimitsWithTrachea(S4R, 1, 0, -10, m_ApexOfTheChest);
+ m_ListOfSupports["S4L"] = LimitsWithTrachea(S4L, 1, 0, -10, m_ApexOfTheChest);
}
//--------------------------------------------------------------------
-//--------------------------------------------------------------------
-template <class ImageType>
-void
-clitk::ExtractLymphStationsFilter<ImageType>::
-Support_S5()
-{
- StartNewStep("[Support] Sup-Inf limits S5 with Aorta and MainPulmonaryArtery");
-
- // Initial S5 support
- MaskImagePointer S5 =
- clitk::Clone<MaskImageType>(this->GetAFDB()->template GetImage<MaskImageType>("Mediastinum", true));
-
- // Sup limits with Aorta
- double sup = FindInferiorBorderOfAorticArch();
-
- // Inf limits with "upper rim of the left main pulmonary artery"
- // For the moment only, it will change.
- MaskImagePointer MainPulmonaryArtery = this->GetAFDB()->template GetImage<MaskImageType>("MainPulmonaryArtery");
- MaskImagePointType p;
- p[0] = p[1] = p[2] = 0.0; // to avoid warning
- clitk::FindExtremaPointInAGivenDirection<MaskImageType>(MainPulmonaryArtery, GetBackgroundValue(), 2, false, p);
- p[2] += MainPulmonaryArtery->GetSpacing()[2];
-
- // Cut Sup/Inf
- S5 = clitk::CropImageAlongOneAxis<MaskImageType>(S5, 2, p[2], sup, true, GetBackgroundValue());
-
- m_ListOfSupports["S5"] = S5;
-}
-//--------------------------------------------------------------------
-
-//--------------------------------------------------------------------
-template <class ImageType>
-void
-clitk::ExtractLymphStationsFilter<ImageType>::
-Support_S6()
-{
- StartNewStep("[Support] Sup-Inf limits S6 with aorta");
-
- // Initial S6 support like S3A
- MaskImagePointer S6 = clitk::Clone<MaskImageType>(m_ListOfSupports["S3A"]);
-
- // Inf Sup limits with Aorta
- double sup = FindSuperiorBorderOfAorticArch();
- double inf = FindInferiorBorderOfAorticArch();
-
- // Cut Sup/Inf
- S6 = clitk::CropImageAlongOneAxis<MaskImageType>(S6, 2, inf, sup, true, GetBackgroundValue());
-
- m_ListOfSupports["S6"] = S6;
-}
-//--------------------------------------------------------------------
-
filter->Update();
} catch(std::runtime_error e) {
std::cout << e.what() << std::endl;
+ return EXIT_FAILURE;
}
return EXIT_SUCCESS;
version "1.0"
purpose "Extract LymphStations with help of TODO"
-option "config" - "Config file" string no
-option "imagetypes" - "Display allowed image types" flag off
-option "verbose" v "Verbose" flag off
-option "verboseStep" - "Verbose each step" flag off
-option "writeStep" w "Write image at each step" flag off
-option "verboseOption" - "Display options values" flag off
-option "verboseWarningOff" - "Do not display warning" flag off
-option "verboseMemory" - "Display memory usage" flag off
+option "config" - "Config file" string no
+option "imagetypes" - "Display allowed image types" flag off
+option "verbose" v "Verbose" flag off
+option "verboseStep" - "Verbose each step" flag off
+option "writeStep" w "Write image at each step" flag off
+option "verboseOption" - "Display options values" flag off
+option "verboseWarningOff" - "Do not display warning" flag off
+option "verboseMemory" - "Display memory usage" flag off
section "I/O"
-option "afdb" a "Input Anatomical Feature DB" string no
-option "input" i "Input filename" string no
-option "output" o "Output lungs mask filename" string no
+option "afdb" a "Input Anatomical Feature DB" string no
+option "afdb_path" p "Input path for image in AFDB" string default="./" no
+option "input" i "Input filename" string no
+option "support_limits" - "Filename to read the support limits" string yes
+option "check_support_limits" - "Display stat on the support limits" flag off
+option "output" o "Output lungs mask filename" string no
section "Options for all stations"
-option "station" - "Force to compute station even if already exist in the DB" string no multiple
-option "nosupport" - "Do not recompute the station supports if already available" flag off
-option "relpos" - "List of filenames for relativeposition operations" string no multiple
+option "force_support" - "Force to compute all supports even if already available" flag off
+option "station" - "Force to compute this station even if already exist in the DB" string no multiple
+option "relpos" - "List of filenames for relativeposition operations" string no multiple
+
+
# section "Options for Station 3A"
# section "Options for Station 2RL"
// clitk
#include "clitkStructuresExtractionFilter.h"
+#include "clitkLabelImageOverlapMeasureFilter.h"
// vtk
#include <vtkPolyData.h>
namespace clitk {
+ class SupportLimitsType {
+ public:
+ std::string station_limit;
+ std::string station;
+ std::string structure_limit;
+ std::string structure;
+ double offset;
+ void Read(istream & is) {
+ is >> station_limit;
+ is >> station;
+ is >> structure_limit;
+ is >> structure;
+ std::string s;
+ is >> s;
+ offset = atof(s.c_str());
+ }
+ };
+
//--------------------------------------------------------------------
/*
Try to extract the LymphStations part of a thorax CT.
double GetFuzzyThreshold(std::string station, std::string tag);
void SetThreshold(std::string station, std::string tag, double value);
double GetThreshold(std::string station, std::string tag);
- itkGetConstMacro(ComputeStationsSupportsFlag, bool);
- itkSetMacro(ComputeStationsSupportsFlag, bool);
- itkBooleanMacro(ComputeStationsSupportsFlag);
+ itkGetConstMacro(ForceSupportsFlag, bool);
+ itkSetMacro(ForceSupportsFlag, bool);
+ itkBooleanMacro(ForceSupportsFlag);
+
+ itkGetConstMacro(CheckSupportFlag, bool);
+ itkSetMacro(CheckSupportFlag, bool);
+ itkBooleanMacro(CheckSupportFlag);
+
+ itkGetConstMacro(SupportLimitsFilename, std::string);
+ itkSetMacro(SupportLimitsFilename, std::string);
protected:
ExtractLymphStationsFilter();
MaskImagePixelType m_BackgroundValue;
MaskImagePixelType m_ForegroundValue;
std::map<std::string, bool> m_ComputeStationMap;
+ std::string m_SupportLimitsFilename;
+ std::vector<SupportLimitsType> m_ListOfSupportLimits;
bool CheckForStation(std::string station);
void Remove_Structures(std::string station, std::string s);
+ void WriteImageSupport(std::string support);
+ void WriteImageStation(std::string station);
+ void ComputeOverlapWithRef(std::string station);
+ void Support_SI_Limit(const std::string station_limit, const std::string station,
+ const std::string structure_limit, const std::string structure,
+ const double offset);
+ void ReadSupportLimits(std::string filename);
// Functions common to several stations
double FindCarina();
// Station's supports
void ExtractStationSupports();
- void Support_SupInf_S1RL();
void Support_LeftRight_S1R_S1L();
- void Support_SupInf_S2R_S2L();
void Support_LeftRight_S2R_S2L();
- void Support_SupInf_S4R_S4L();
void Support_LeftRight_S4R_S4L();
void Support_Post_S1S2S4();
- void Support_S3P();
- void Support_S3A();
- void Support_S5();
- void Support_S6();
MaskImagePointer LimitsWithTrachea(MaskImageType * input,
int extremaDirection, int lineDirection,
// Station 4RL
void ExtractStation_4RL_SetDefaultValues();
- void ExtractStation_4RL();
+ void ExtractStation_4L();
+ void ExtractStation_4R();
void ExtractStation_S4L_S5_Limits_Aorta_LeftPulmonaryArtery(int KeepPoint);
// Station 5
void ExtractStation_7_Posterior_Limits();
void ExtractStation_7_Remove_Structures();
bool m_S7_UseMostInferiorPartOnlyFlag;
- bool m_ComputeStationsSupportsFlag;
+ bool m_ForceSupportsFlag;
+ bool m_CheckSupportFlag;
MaskImagePointer m_Working_Trachea;
MaskImagePointer m_LeftBronchus;
MaskImagePointer m_RightBronchus;
+++ /dev/null
-/*=========================================================================
- Program: vv http://www.creatis.insa-lyon.fr/rio/vv
-
- Authors belong to:
- - University of LYON http://www.universite-lyon.fr/
- - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr
- - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
-
- This software is distributed WITHOUT ANY WARRANTY; without even
- the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
- PURPOSE. See the copyright notices for more information.
-
- It is distributed under dual licence
-
- - BSD See included LICENSE.txt file
- - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
- ======================================================================-====*/
-
-#ifndef CLITKEXTRACTLYMPHSTATIONSFILTER_H
-#define CLITKEXTRACTLYMPHSTATIONSFILTER_H
-
-// clitk
-#include "clitkFilterBase.h"
-#include "clitkFilterWithAnatomicalFeatureDatabaseManagement.h"
-
-namespace clitk {
-
- //--------------------------------------------------------------------
- /*
- Try to extract the LymphStations part of a thorax CT.
- Inputs :
- - Patient label image
- - Lungs label image
- - Bones label image
- */
- //--------------------------------------------------------------------
-
- template <class TImageType>
- class ITK_EXPORT ExtractLymphStationsFilter:
- public virtual clitk::FilterBase,
- public clitk::FilterWithAnatomicalFeatureDatabaseManagement,
- public itk::ImageToImageFilter<TImageType, itk::Image<uchar, 3> >
- {
-
- public:
- /** Standard class typedefs. */
- typedef itk::ImageToImageFilter<TImageType, itk::Image<uchar, 3> > Superclass;
- typedef ExtractLymphStationsFilter Self;
- typedef itk::SmartPointer<Self> Pointer;
- typedef itk::SmartPointer<const Self> ConstPointer;
-
- /** Method for creation through the object factory. */
- itkNewMacro(Self);
-
- /** Run-time type information (and related methods). */
- itkTypeMacro(ExtractLymphStationsFilter, InPlaceImageFilter);
-
- /** Some convenient typedefs. */
- typedef TImageType ImageType;
- typedef typename ImageType::ConstPointer ImageConstPointer;
- typedef typename ImageType::Pointer ImagePointer;
- typedef typename ImageType::RegionType ImageRegionType;
- typedef typename ImageType::PixelType ImagePixelType;
- typedef typename ImageType::SizeType ImageSizeType;
- typedef typename ImageType::IndexType ImageIndexType;
- typedef typename ImageType::PointType ImagePointType;
-
- typedef uchar MaskImagePixelType;
- typedef itk::Image<MaskImagePixelType, 3> MaskImageType;
- typedef typename MaskImageType::Pointer MaskImagePointer;
-
- /** Connect inputs */
- // void SetInput(const TImageType * image);
-
- /** ImageDimension constants */
- itkStaticConstMacro(ImageDimension, unsigned int, TImageType::ImageDimension);
- FILTERBASE_INIT;
-
- itkGetConstMacro(BackgroundValue, ImagePixelType);
- itkGetConstMacro(ForegroundValue, ImagePixelType);
- itkSetMacro(BackgroundValue, ImagePixelType);
- itkSetMacro(ForegroundValue, ImagePixelType);
-
- itkSetMacro(FuzzyThreshold, double);
- itkGetConstMacro(FuzzyThreshold, double);
-
- protected:
- ExtractLymphStationsFilter();
- virtual ~ExtractLymphStationsFilter() {}
-
- virtual void GenerateOutputInformation();
- virtual void GenerateInputRequestedRegion();
- virtual void GenerateData();
-
- MaskImagePointer m_mediastinum;
- MaskImagePointer m_trachea;
- MaskImagePointer m_working_mediastinum;
- MaskImagePointer m_working_trachea;
- MaskImagePointer m_output;
-
- ImagePixelType m_BackgroundValue;
- ImagePixelType m_ForegroundValue;
- double m_CarinaZPositionInMM;
- bool m_CarinaZPositionInMMIsSet;
- double m_MiddleLobeBronchusZPositionInMM;
-
- double m_IntermediateSpacing;
- double m_FuzzyThreshold;
-
- private:
- ExtractLymphStationsFilter(const Self&); //purposely not implemented
- void operator=(const Self&); //purposely not implemented
-
- }; // end class
- //--------------------------------------------------------------------
-
-} // end namespace clitk
-//--------------------------------------------------------------------
-
-#ifndef ITK_MANUAL_INSTANTIATION
-#include "clitkExtractLymphStationsFilter.txx"
-#endif
-
-#endif
+++ /dev/null
-/*=========================================================================
- Program: vv http://www.creatis.insa-lyon.fr/rio/vv
-
- Authors belong to:
- - University of LYON http://www.universite-lyon.fr/
- - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr
- - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
-
- This software is distributed WITHOUT ANY WARRANTY; without even
- the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
- PURPOSE. See the copyright notices for more information.
-
- It is distributed under dual licence
-
- - BSD See included LICENSE.txt file
- - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
- ======================================================================-====*/
-
-#ifndef CLITKEXTRACTLYMPHSTATIONSFILTER_TXX
-#define CLITKEXTRACTLYMPHSTATIONSFILTER_TXX
-
-// clitk
-#include "clitkCommon.h"
-#include "clitkExtractLymphStationsFilter.h"
-#include "clitkAddRelativePositionConstraintToLabelImageFilter.h"
-#include "clitkSegmentationUtils.h"
-#include "clitkAutoCropFilter.h"
-#include "clitkSegmentationUtils.h"
-#include "clitkSliceBySliceRelativePositionFilter.h"
-
-// itk
-#include <itkStatisticsLabelMapFilter.h>
-#include <itkLabelImageToStatisticsLabelMapFilter.h>
-#include <itkRegionOfInterestImageFilter.h>
-#include <itkBinaryThresholdImageFilter.h>
-#include <itkImageSliceConstIteratorWithIndex.h>
-#include <itkImageSliceIteratorWithIndex.h>
-#include <itkBinaryThinningImageFilter.h>
-
-// itk ENST
-#include "RelativePositionPropImageFilter.h"
-
-//--------------------------------------------------------------------
-template <class TImageType>
-clitk::ExtractLymphStationsFilter<TImageType>::
-ExtractLymphStationsFilter():
- clitk::FilterBase(),
- clitk::FilterWithAnatomicalFeatureDatabaseManagement(),
- itk::ImageToImageFilter<TImageType, MaskImageType>()
-{
- this->SetNumberOfRequiredInputs(1);
- SetBackgroundValue(0);
- SetForegroundValue(1);
- SetFuzzyThreshold(0.5);
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class TImageType>
-void
-clitk::ExtractLymphStationsFilter<TImageType>::
-GenerateOutputInformation() {
- // Superclass::GenerateOutputInformation();
-
- // Get input
- LoadAFDB();
- m_mediastinum = GetAFDB()->template GetImage <MaskImageType>("mediastinum");
- m_trachea = GetAFDB()->template GetImage <MaskImageType>("trachea");
- ImagePointType carina;
- GetAFDB()->GetPoint3D("carina", carina);
- m_CarinaZPositionInMM = carina[2];
- DD(m_CarinaZPositionInMM);
- ImagePointType mlb;
- GetAFDB()->GetPoint3D("middleLobeBronchus", mlb);
- m_MiddleLobeBronchusZPositionInMM = mlb[2];
- DD(m_MiddleLobeBronchusZPositionInMM);
-
- // ----------------------------------------------------------------
- // ----------------------------------------------------------------
- // Superior limit = carina
- // Inferior limit = origin right middle lobe bronchus
- StartNewStep("Inf/Sup mediastinum limits with carina/bronchus");
- m_working_mediastinum =
- clitk::CropImageAlongOneAxis<MaskImageType>(m_mediastinum, 2,
- m_MiddleLobeBronchusZPositionInMM,
- m_CarinaZPositionInMM, true,
- GetBackgroundValue());
- StopCurrentStep<MaskImageType>(m_working_mediastinum);
- m_output = m_working_mediastinum;
-
- // ----------------------------------------------------------------
- // ----------------------------------------------------------------
- // Superior limit = carina
- // Inferior limit = origin right middle lobe bronchus
- StartNewStep("Inf/Sup trachea limits with carina/bronchus");
- m_working_trachea =
- clitk::CropImageAlongOneAxis<MaskImageType>(m_trachea, 2,
- m_MiddleLobeBronchusZPositionInMM,
- m_CarinaZPositionInMM, true,
- GetBackgroundValue());
- StopCurrentStep<MaskImageType>(m_working_trachea);
-
- // ----------------------------------------------------------------
- // ----------------------------------------------------------------
- // Separate trachea in two CCL
- StartNewStep("Separate trachea under carina");
-
- // Labelize and consider two main labels
- m_working_trachea = Labelize<MaskImageType>(m_working_trachea, 0, true, 1);
-
- // Carina position must at the first slice that separate the two main bronchus (not superiorly)
- // Check that upper slice is composed of at least two labels
- typedef itk::ImageSliceIteratorWithIndex<MaskImageType> SliceIteratorType;
- SliceIteratorType iter(m_working_trachea, m_working_trachea->GetLargestPossibleRegion());
- iter.SetFirstDirection(0);
- iter.SetSecondDirection(1);
- iter.GoToReverseBegin(); // Start from the end (because image is IS not SI)
- int maxLabel=0;
- while (!iter.IsAtReverseEndOfSlice()) {
- while (!iter.IsAtReverseEndOfLine()) {
- // DD(iter.GetIndex());
- if (iter.Get() > maxLabel) maxLabel = iter.Get();
- --iter;
- }
- iter.PreviousLine();
- }
- if (maxLabel < 2) {
- clitkExceptionMacro("First slice form Carina does not seems to seperate the two main bronchus. Abort");
- }
-
- // Compute centroid of both parts to identify the left from the right bronchus
- typedef long LabelType;
- static const unsigned int Dim = ImageType::ImageDimension;
- typedef itk::ShapeLabelObject< LabelType, Dim > LabelObjectType;
- typedef itk::LabelMap< LabelObjectType > LabelMapType;
- typedef itk::LabelImageToLabelMapFilter<MaskImageType, LabelMapType> ImageToMapFilterType;
- typename ImageToMapFilterType::Pointer imageToLabelFilter = ImageToMapFilterType::New();
- typedef itk::ShapeLabelMapFilter<LabelMapType, MaskImageType> ShapeFilterType;
- typename ShapeFilterType::Pointer statFilter = ShapeFilterType::New();
- imageToLabelFilter->SetBackgroundValue(GetBackgroundValue());
- imageToLabelFilter->SetInput(m_working_trachea);
- statFilter->SetInput(imageToLabelFilter->GetOutput());
- statFilter->Update();
- typename LabelMapType::Pointer labelMap = statFilter->GetOutput();
-
- ImagePointType C1 = labelMap->GetLabelObject(1)->GetCentroid();
- ImagePointType C2 = labelMap->GetLabelObject(2)->GetCentroid();
-
- ImagePixelType leftLabel;
- ImagePixelType rightLabel;
- if (C1[0] < C2[0]) { leftLabel = 1; rightLabel = 2; }
- else { leftLabel = 2; rightLabel = 1; }
- DD(leftLabel);
- DD(rightLabel);
-
- StopCurrentStep<MaskImageType>(m_working_trachea);
-
- //-----------------------------------------------------
- StartNewStep("Left limits with bronchus (slice by slice)");
- // Select LeftLabel (set one label to Backgroundvalue)
- MaskImagePointer leftBronchus =
- SetBackground<MaskImageType, MaskImageType>(m_working_trachea, m_working_trachea,
- rightLabel, GetBackgroundValue());
- MaskImagePointer rightBronchus =
- SetBackground<MaskImageType, MaskImageType>(m_working_trachea, m_working_trachea,
- leftLabel, GetBackgroundValue());
- writeImage<MaskImageType>(leftBronchus, "left.mhd");
- writeImage<MaskImageType>(rightBronchus, "right.mhd");
-
- m_working_mediastinum =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_working_mediastinum,
- leftBronchus, 2,
- GetFuzzyThreshold(), "RightTo",
- true, 4);
- m_working_mediastinum =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_working_mediastinum,
- rightBronchus,
- 2, GetFuzzyThreshold(), "LeftTo",
- true, 4);
- m_working_mediastinum =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_working_mediastinum,
- leftBronchus,
- 2, GetFuzzyThreshold(), "AntTo",
- true, 4, true); // NOT
- m_working_mediastinum =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_working_mediastinum,
- rightBronchus,
- 2, GetFuzzyThreshold(), "AntTo",
- true, 4, true);
- m_working_mediastinum =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_working_mediastinum,
- leftBronchus,
- 2, GetFuzzyThreshold(), "PostTo",
- true, 4, true);
- m_working_mediastinum =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_working_mediastinum,
- rightBronchus,
- 2, GetFuzzyThreshold(), "PostTo",
- true, 4, true);
- m_output = m_working_mediastinum;
- StopCurrentStep<MaskImageType>(m_output);
-
- //-----------------------------------------------------
- //-----------------------------------------------------
- StartNewStep("Posterior limits");
-
- // Left Bronchus slices
- typedef clitk::ExtractSliceFilter<MaskImageType> ExtractSliceFilterType;
- typedef typename ExtractSliceFilterType::SliceType SliceType;
- typename ExtractSliceFilterType::Pointer
- extractSliceFilter = ExtractSliceFilterType::New();
- extractSliceFilter->SetInput(leftBronchus);
- extractSliceFilter->SetDirection(2);
- extractSliceFilter->Update();
- std::vector<typename SliceType::Pointer> leftBronchusSlices;
- extractSliceFilter->GetOutputSlices(leftBronchusSlices);
-
- // Right Bronchus slices
- extractSliceFilter = ExtractSliceFilterType::New();
- extractSliceFilter->SetInput(rightBronchus);
- extractSliceFilter->SetDirection(2);
- extractSliceFilter->Update();
- std::vector<typename SliceType::Pointer> rightBronchusSlices ;
- extractSliceFilter->GetOutputSlices(rightBronchusSlices);
-
- assert(leftBronchusSlices.size() == rightBronchusSlices.size());
-
- std::vector<MaskImageType::PointType> leftPoints;
- std::vector<MaskImageType::PointType> rightPoints;
- for(uint i=0; i<leftBronchusSlices.size(); i++) {
- // Keep main CCL
- leftBronchusSlices[i] = Labelize<SliceType>(leftBronchusSlices[i], 0, true, 10);
- leftBronchusSlices[i] = KeepLabels<SliceType>(leftBronchusSlices[i],
- GetBackgroundValue(),
- GetForegroundValue(), 1, 1, true);
- rightBronchusSlices[i] = Labelize<SliceType>(rightBronchusSlices[i], 0, true, 10);
- rightBronchusSlices[i] = KeepLabels<SliceType>(rightBronchusSlices[i],
- GetBackgroundValue(),
- GetForegroundValue(), 1, 1, true);
- double distance_max_from_center_point = 15;
-
- // ------- Find point in left Bronchus -------
- // find right most point in left = rightMost
- SliceType::PointType a;
- SliceType::PointType rightMost =
- clitk::FindExtremaPointInAGivenDirection<SliceType>(leftBronchusSlices[i],
- GetBackgroundValue(),
- 0, false, a, 0);
- // find post most point in left, not far away from rightMost
- SliceType::PointType p =
- clitk::FindExtremaPointInAGivenDirection<SliceType>(leftBronchusSlices[i],
- GetBackgroundValue(),
- 1, false, rightMost,
- distance_max_from_center_point);
- MaskImageType::PointType pp;
- pp[0] = p[0]; pp[1] = p[1];
- pp[2] = i*leftBronchus->GetSpacing()[2] + leftBronchus->GetOrigin()[2];
- leftPoints.push_back(pp);
-
- // ------- Find point in right Bronchus -------
- // find left most point in right = leftMost
- SliceType::PointType leftMost =
- clitk::FindExtremaPointInAGivenDirection<SliceType>(rightBronchusSlices[i],
- GetBackgroundValue(),
- 0, true, a, 0);
- // find post most point in left, not far away from leftMost
- p = clitk::FindExtremaPointInAGivenDirection<SliceType>(rightBronchusSlices[i],
- GetBackgroundValue(),
- 1, false, leftMost,
- distance_max_from_center_point);
- pp[0] = p[0]; pp[1] = p[1];
- pp[2] = i*rightBronchus->GetSpacing()[2] + rightBronchus->GetOrigin()[2];
- rightPoints.push_back(pp);
- }
-
- // DEBUG
- std::ofstream osl;
- openFileForWriting(osl, "left.txt");
- osl << "LANDMARKS1" << std::endl;
- std::ofstream osr;
- openFileForWriting(osr, "right.txt");
- osr << "LANDMARKS1" << std::endl;
- for(uint i=0; i<leftBronchusSlices.size(); i++) {
- osl << i << " " << leftPoints[i][0] << " " << leftPoints[i][1]
- << " " << leftPoints[i][2] << " 0 0 " << std::endl;
- osr << i << " " << rightPoints[i][0] << " " << rightPoints[i][1]
- << " " << rightPoints[i][2] << " 0 0 " << std::endl;
- }
- osl.close();
- osr.close();
-
- // Now uses these points to limit, slice by slice
- // http://www.gamedev.net/community/forums/topic.asp?topic_id=542870
- /*
- Assuming the points are (Ax,Ay) (Bx,By) and (Cx,Cy), you need to compute:
- (Bx - Ax) * (Cy - Ay) - (By - Ay) * (Cx - Ax)
- This will equal zero if the point C is on the line formed by points A and B, and will have a different sign depending on the side. Which side this is depends on the orientation of your (x,y) coordinates, but you can plug test values for A,B and C into this formula to determine whether negative values are to the left or to the right.
-
- => to accelerate, start with formula, when change sign -> stop and fill
- */
- // typedef itk::ImageSliceIteratorWithIndex<MaskImageType> SliceIteratorType;
- iter = SliceIteratorType(m_working_mediastinum, m_working_mediastinum->GetLargestPossibleRegion());
- iter.SetFirstDirection(0);
- iter.SetSecondDirection(1);
- iter.GoToBegin();
- int i=0;
- MaskImageType::PointType A;
- MaskImageType::PointType B;
- MaskImageType::PointType C;
- while (!iter.IsAtEnd()) {
- A = leftPoints[i];
- B = rightPoints[i];
- C = A;
- C[1] -= 10; // I know I must keep this point
- double s = (B[0] - A[0]) * (C[1] - A[1]) - (B[1] - A[1]) * (C[0] - A[0]);
- bool isPositive = s<0;
- while (!iter.IsAtEndOfSlice()) {
- while (!iter.IsAtEndOfLine()) {
- // Very slow, I know ... but image should be very small
- m_working_mediastinum->TransformIndexToPhysicalPoint(iter.GetIndex(), C);
- double s = (B[0] - A[0]) * (C[1] - A[1]) - (B[1] - A[1]) * (C[0] - A[0]);
- if (s == 0) iter.Set(2);
- if (isPositive) {
- if (s > 0) iter.Set(GetBackgroundValue());
- }
- else {
- if (s < 0) iter.Set(GetBackgroundValue());
- }
- ++iter;
- }
- iter.NextLine();
- }
- iter.NextSlice();
- ++i;
- }
-
- //-----------------------------------------------------
- //-----------------------------------------------------
- // StartNewStep("Anterior limits");
-
-
- // MISSING FROM NOW
-
- // Station 4R, Station 4L, the right pulmonary artery, and/or the
- // left superior pulmonary vein
-
-
- //-----------------------------------------------------
- //-----------------------------------------------------
- // ALSO SUBSTRACT ARTERY/VEIN (initially in the support)
-
-
- // Set output image information (required)
- MaskImagePointer outputImage = this->GetOutput(0);
- outputImage->SetRegions(m_working_mediastinum->GetLargestPossibleRegion());
- outputImage->SetOrigin(m_working_mediastinum->GetOrigin());
- outputImage->SetRequestedRegion(m_working_mediastinum->GetLargestPossibleRegion());
-
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class TImageType>
-void
-clitk::ExtractLymphStationsFilter<TImageType>::
-GenerateInputRequestedRegion() {
- DD("GenerateInputRequestedRegion");
- // Call default
- // Superclass::GenerateInputRequestedRegion();
- // Following needed because output region can be greater than input (trachea)
- //ImagePointer mediastinum = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(0));
- //ImagePointer trachea = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(1));
- DD("set reg");
- m_mediastinum->SetRequestedRegion(m_mediastinum->GetLargestPossibleRegion());
- m_trachea->SetRequestedRegion(m_trachea->GetLargestPossibleRegion());
- DD("end");
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class TImageType>
-void
-clitk::ExtractLymphStationsFilter<TImageType>::
-GenerateData() {
- DD("GenerateData");
- // Final Step -> graft output (if SetNthOutput => redo)
- this->GraftOutput(m_output);
-}
-//--------------------------------------------------------------------
-
-
-#endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX
this->SetNumberOfRequiredInputs(1);
SetBackgroundValue(0);
SetForegroundValue(1);
- ComputeStationsSupportsFlagOn();
+ ForceSupportsFlagOn();
+ SetSupportLimitsFilename("none");
+ CheckSupportFlagOff();
// Default values
ExtractStation_3P_SetDefaultValues();
m_Mediastinum = this->GetAFDB()->template GetImage <MaskImageType>("Mediastinum");
// Clean some computer landmarks to force the recomputation
+ // FIXME -> to put elsewhere ?
this->GetAFDB()->RemoveTag("AntPostVesselsSeparation");
- // Global supports for stations
+ // Must I compute the supports ?
bool supportsExist = true;
- try {
- m_ListOfSupports["S1R"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S1R");
- m_ListOfSupports["S1L"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S1L");
- m_ListOfSupports["S2R"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S2R");
- m_ListOfSupports["S2L"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S2L");
- m_ListOfSupports["S4R"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S4R");
- m_ListOfSupports["S4L"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S4L");
-
- m_ListOfSupports["S3A"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S3A");
- m_ListOfSupports["S3P"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S3P");
- m_ListOfSupports["S5"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S5");
- m_ListOfSupports["S6"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S6");
- m_ListOfSupports["S7"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S7");
- m_ListOfSupports["S8"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S8");
- m_ListOfSupports["S9"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S9");
- m_ListOfSupports["S10"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S10");
- m_ListOfSupports["S11"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S11");
- } catch(clitk::ExceptionObject o) {
- supportsExist = false;
+ if (!GetForceSupportsFlag()) {
+ try {
+ m_ListOfSupports["S1R"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S1R");
+ m_ListOfSupports["S1L"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S1L");
+ m_ListOfSupports["S2R"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S2R");
+ m_ListOfSupports["S2L"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S2L");
+ m_ListOfSupports["S4R"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S4R");
+ m_ListOfSupports["S4L"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S4L");
+
+ m_ListOfSupports["S3A"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S3A");
+ m_ListOfSupports["S3P"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S3P");
+ m_ListOfSupports["S5"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S5");
+ m_ListOfSupports["S6"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S6");
+ m_ListOfSupports["S7"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S7");
+ m_ListOfSupports["S8"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S8");
+ m_ListOfSupports["S9"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S9");
+ m_ListOfSupports["S10"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S10");
+ m_ListOfSupports["S11"] = this->GetAFDB()->template GetImage<MaskImageType>("Support_S11");
+ } catch(clitk::ExceptionObject o) {
+ supportsExist = false;
+ }
}
- if (!supportsExist || GetComputeStationsSupportsFlag()) {
+ if (!supportsExist || GetForceSupportsFlag()) {
this->StartNewStep("Supports for stations");
this->StartSubStep();
+
+ // FIXME : why should I remove theses tags ???
this->GetAFDB()->RemoveTag("CarinaZ");
this->GetAFDB()->RemoveTag("ApexOfTheChestZ");
this->GetAFDB()->RemoveTag("ApexOfTheChest");
ExtractStation_2RL();
ExtractStation_3P();
ExtractStation_3A();
- ExtractStation_4RL();
+ ExtractStation_4R();
+ ExtractStation_4L();
ExtractStation_5();
ExtractStation_6();
- // ---------- TODO -----------------------
+ // ---------- todo -----------------------
// Extract Station8
// ExtractStation_8();
clitk::ExtractLymphStationsFilter<TImageType>::
GenerateData() {
// Final Step -> graft output (if SetNthOutput => redo)
- this->GraftOutput(m_ListOfStations["8"]);
+ // this->GraftOutput(m_ListOfStations["8"]);
}
//--------------------------------------------------------------------
std::string s = "Station"+station;
+ // Define the starting support
+ // if (GetComputeStation(station)) {
+ // std::cout << "Station " << station << " already exists, but re-computation forced." << std::endl;
+ // }
+ if (GetComputeStation(station)) {
+ m_Working_Support = m_Mediastinum = this->GetAFDB()->template GetImage<MaskImageType>("Mediastinum", true);
+ return true;
+ }
+ else return false;
+ // else {
+ // std::cout << "Station " << station << " found. I used it" << std::endl;
+ // return false;
+ // }
+
// Check if station already exist in DB
+
+ // FIXME -> do nothing if not on the command line. Is it what I want ?
bool found = false;
if (this->GetAFDB()->TagExist(s)) {
m_ListOfStations[station] = this->GetAFDB()->template GetImage<MaskImageType>(s);
found = true;
}
- // Define the starting support
- if (found && GetComputeStation(station)) {
- std::cout << "Station " << station << " already exists, but re-computation forced." << std::endl;
- }
- if (!found || GetComputeStation(station)) {
- m_Working_Support = m_Mediastinum = this->GetAFDB()->template GetImage<MaskImageType>("Mediastinum", true);
- return true;
- }
- else {
- std::cout << "Station " << station << " found. I used it" << std::endl;
- return false;
- }
}
//--------------------------------------------------------------------
z = this->GetAFDB()->GetDouble("CarinaZ");
}
catch(clitk::ExceptionObject e) {
- DD("FindCarinaSlicePosition");
+ //DD("FindCarinaSlicePosition");
// Get Carina
MaskImagePointer Carina = this->GetAFDB()->template GetImage<MaskImageType>("Carina");
z = this->GetAFDB()->GetDouble("ApexOfTheChestZ");
}
catch(clitk::ExceptionObject e) {
- DD("FindApexOfTheChestPosition");
+
+ /*
+ //DD("FindApexOfTheChestPosition");
MaskImagePointer Lungs = this->GetAFDB()->template GetImage<MaskImageType>("Lungs");
MaskImagePointType p;
p[0] = p[1] = p[2] = 0.0; // to avoid warning
this->GetAFDB()->SetDouble("ApexOfTheChestZ", p[2]);
this->WriteAFDB();
z = p[2];
+ */
+
+ // the superior border becomes the more inferior of the two apices
+ MaskImagePointer RightLung = this->GetAFDB()->template GetImage<MaskImageType>("RightLung");
+ MaskImagePointer LeftLung = this->GetAFDB()->template GetImage<MaskImageType>("LeftLung");
+ MaskImagePointType pr;
+ MaskImagePointType pl;
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(RightLung, GetBackgroundValue(), 2, false, pr);
+ clitk::FindExtremaPointInAGivenDirection<MaskImageType>(LeftLung, GetBackgroundValue(), 2, false, pl);
+ // We dont need Lungs structure from now
+ this->GetAFDB()->template ReleaseImage<MaskImageType>("RightLung");
+ this->GetAFDB()->template ReleaseImage<MaskImageType>("LeftLung");
+ // Put inside the AFDB
+ if (pr[2] < pl[2]) z = pr[2];
+ else z = pl[2];
+ this->GetAFDB()->SetDouble("ApexOfTheChestZ", z);
+ this->WriteAFDB();
}
return z;
}
z = this->GetAFDB()->GetDouble("SuperiorBorderOfAorticArchZ");
}
catch(clitk::ExceptionObject e) {
- DD("FindSuperiorBorderOfAorticArch");
+ // DD("FindSuperiorBorderOfAorticArch");
MaskImagePointer Aorta = this->GetAFDB()->template GetImage<MaskImageType>("Aorta");
MaskImagePointType p;
p[0] = p[1] = p[2] = 0.0; // to avoid warning
z = this->GetAFDB()->GetDouble("InferiorBorderOfAorticArchZ");
}
catch(clitk::ExceptionObject e) {
- DD("FindInferiorBorderOfAorticArch");
+ //DD("FindInferiorBorderOfAorticArch");
MaskImagePointer Aorta = this->GetAFDB()->template GetImage<MaskImageType>("Aorta");
std::vector<MaskSlicePointer> slices;
clitk::ExtractSlices<MaskImageType>(Aorta, 2, slices);
}
//--------------------------------------------------------------------
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractLymphStationsFilter<ImageType>::
+WriteImageSupport(std::string support)
+{
+ writeImage<MaskImageType>(m_ListOfSupports[support], this->GetAFDBPath()+"/"+"seg/Support_"+support+".mha");
+ this->GetAFDB()->SetImageFilename("Support_"+support, "seg/Support_"+support+".mha");
+}
+//--------------------------------------------------------------------
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractLymphStationsFilter<ImageType>::
+WriteImageStation(std::string station)
+{
+ writeImage<MaskImageType>(m_ListOfStations[station], GetAFDB()->GetPath()+"/seg/Station"+station+".mha");
+ GetAFDB()->SetImageFilename("Station"+station, "seg/Station"+station+".mha");
+ WriteAFDB();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractLymphStationsFilter<ImageType>::
+ComputeOverlapWithRef(std::string station)
+{
+ if (GetComputeStation(station)) {
+ MaskImagePointer ref = this->GetAFDB()->template GetImage <MaskImageType>("Station"+station+"_Ref");
+ typedef clitk::LabelImageOverlapMeasureFilter<MaskImageType> FilterType;
+ typename FilterType::Pointer filter = FilterType::New();
+ filter->SetInput(0, m_ListOfStations[station]);
+ filter->SetInput(1, ref);
+ filter->Update();
+ }
+}
+//--------------------------------------------------------------------
+
+//--------------------------------------------------------------------
+template <class ImageType>
+void
+clitk::ExtractLymphStationsFilter<ImageType>::
+ReadSupportLimits(std::string filename)
+{
+ m_ListOfSupportLimits.clear();
+ ifstream is;
+ openFileForReading(is, filename);
+ while (is) {
+ skipComment(is);
+ SupportLimitsType s;
+ s.Read(is);
+ if (is) m_ListOfSupportLimits.push_back(s);
+ }
+}
+//--------------------------------------------------------------------
+
#endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX
+
+++ /dev/null
-/*=========================================================================
- Program: vv http://www.creatis.insa-lyon.fr/rio/vv
-
- Authors belong to:
- - University of LYON http://www.universite-lyon.fr/
- - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr
- - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
-
- This software is distributed WITHOUT ANY WARRANTY; without even
- the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
- PURPOSE. See the copyright notices for more information.
-
- It is distributed under dual licence
-
- - BSD See included LICENSE.txt file
- - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
- ======================================================================-====*/
-
-#ifndef CLITKEXTRACTLYMPHSTATIONSFILTER_TXX
-#define CLITKEXTRACTLYMPHSTATIONSFILTER_TXX
-
-// clitk
-#include "clitkCommon.h"
-#include "clitkExtractLymphStationsFilter.h"
-#include "clitkAddRelativePositionConstraintToLabelImageFilter.h"
-#include "clitkSegmentationUtils.h"
-#include "clitkAutoCropFilter.h"
-#include "clitkSegmentationUtils.h"
-#include "clitkSliceBySliceRelativePositionFilter.h"
-
-// itk
-#include <itkStatisticsLabelMapFilter.h>
-#include <itkLabelImageToStatisticsLabelMapFilter.h>
-#include <itkRegionOfInterestImageFilter.h>
-#include <itkBinaryThresholdImageFilter.h>
-#include <itkImageSliceConstIteratorWithIndex.h>
-#include <itkBinaryThinningImageFilter.h>
-
-// itk ENST
-#include "RelativePositionPropImageFilter.h"
-
-//--------------------------------------------------------------------
-template <class TImageType>
-clitk::ExtractLymphStationsFilter<TImageType>::
-ExtractLymphStationsFilter():
- clitk::FilterBase(),
- clitk::FilterWithAnatomicalFeatureDatabaseManagement(),
- itk::ImageToImageFilter<TImageType, MaskImageType>()
-{
- this->SetNumberOfRequiredInputs(1);
- SetBackgroundValue(0);
- SetForegroundValue(1);
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class TImageType>
-void
-clitk::ExtractLymphStationsFilter<TImageType>::
-GenerateOutputInformation() {
- // Superclass::GenerateOutputInformation();
-
- // Get input
- LoadAFDB();
- m_mediastinum = GetAFDB()->template GetImage <MaskImageType>("mediastinum");
- m_trachea = GetAFDB()->template GetImage <MaskImageType>("trachea");
- ImagePointType carina;
- GetAFDB()->GetPoint3D("carina", carina);
- DD(carina);
- m_CarinaZPositionInMM = carina[2];
- DD(m_CarinaZPositionInMM);
- ImagePointType mlb;
- GetAFDB()->GetPoint3D("middleLobeBronchus", mlb);
- m_MiddleLobeBronchusZPositionInMM = mlb[2];
- DD(m_MiddleLobeBronchusZPositionInMM);
-
- // ----------------------------------------------------------------
- // ----------------------------------------------------------------
- // Superior limit = carina
- // Inferior limit = origin right middle lobe bronchus
- StartNewStep("Inf/Sup mediastinum limits with carina/bronchus");
- ImageRegionType region = m_mediastinum->GetLargestPossibleRegion(); DD(region);
- ImageSizeType size = region.GetSize();
- ImageIndexType index = region.GetIndex();
- DD(m_CarinaZPositionInMM);
- DD(m_MiddleLobeBronchusZPositionInMM);
- index[2] = floor((m_MiddleLobeBronchusZPositionInMM - m_mediastinum->GetOrigin()[2]) / m_mediastinum->GetSpacing()[2]);
- size[2] = 1+ceil((m_CarinaZPositionInMM-m_MiddleLobeBronchusZPositionInMM) / m_mediastinum->GetSpacing()[2]); // +1 to
- region.SetSize(size);
- region.SetIndex(index);
- DD(region);
- typedef itk::RegionOfInterestImageFilter<MaskImageType, MaskImageType> CropFilterType;
- typename CropFilterType::Pointer cropFilter = CropFilterType::New();
- cropFilter->SetInput(m_mediastinum);
- cropFilter->SetRegionOfInterest(region);
- cropFilter->Update();
- m_working_image = cropFilter->GetOutput();
- // Auto Crop (because following rel pos is faster)
- m_working_image = clitk::AutoCrop<MaskImageType>(m_working_image, 0);
- StopCurrentStep<MaskImageType>(m_working_image);
- m_output = m_working_image;
-
- // ----------------------------------------------------------------
- // ----------------------------------------------------------------
- // Separate trachea in two CCL
- StartNewStep("Separate trachea under carina");
- // DD(region);
- ImageRegionType trachea_region = m_trachea->GetLargestPossibleRegion();
- for(int i=0; i<3; i++) {
- index[i] = floor(((index[i]*m_mediastinum->GetSpacing()[i])+m_mediastinum->GetOrigin()[i]
- -m_trachea->GetOrigin()[i])/m_trachea->GetSpacing()[i]);
- // DD(index[i]);
- size[i] = ceil((size[i]*m_mediastinum->GetSpacing()[i])/m_trachea->GetSpacing()[i]);
- // DD(size[i]);
- if (index[i] < 0) {
- size[i] += index[i];
- index[i] = 0;
- }
- if (size[i]+index[i] > (trachea_region.GetSize()[i] + trachea_region.GetIndex()[i])) {
- size[i] = trachea_region.GetSize()[i] + trachea_region.GetIndex()[i] - index[i];
- }
- }
- // DD(index);
- // DD(size);
- region.SetIndex(index);
- region.SetSize(size);
- // typedef itk::RegionOfInterestImageFilter<MaskImageType, MaskImageType> CropFilterType;
- // typename CropFilterType::Pointer
- cropFilter = CropFilterType::New();
- // m_trachea.Print(std::cout);
- cropFilter->SetInput(m_trachea);
- cropFilter->SetRegionOfInterest(region);
- cropFilter->Update();
- m_working_trachea = cropFilter->GetOutput();
-
- // Labelize and consider two main labels
- m_working_trachea = Labelize<MaskImageType>(m_working_trachea, 0, true, 1);
-
- // Detect wich label is at Left
- typedef itk::ImageSliceConstIteratorWithIndex<MaskImageType> SliceIteratorType;
- SliceIteratorType iter(m_working_trachea, m_working_trachea->GetLargestPossibleRegion());
- iter.SetFirstDirection(0);
- iter.SetSecondDirection(1);
- iter.GoToBegin();
- bool stop = false;
- ImagePixelType leftLabel;
- ImagePixelType rightLabel;
- while (!stop) {
- if (iter.Get() != GetBackgroundValue()) {
- // DD(iter.GetIndex());
- // DD((int)iter.Get());
- leftLabel = iter.Get();
- stop = true;
- }
- ++iter;
- }
- if (leftLabel == 1) rightLabel = 2;
- else rightLabel = 1;
- DD((int)leftLabel);
- DD((int)rightLabel);
-
- // End step
- StopCurrentStep<MaskImageType>(m_working_trachea);
-
- //-----------------------------------------------------
- /* DD("TEST Skeleton");
- typedef itk::BinaryThinningImageFilter<MaskImageType, MaskImageType> SkeletonFilterType;
- typename SkeletonFilterType::Pointer skeletonFilter = SkeletonFilterType::New();
- skeletonFilter->SetInput(m_working_trachea);
- skeletonFilter->Update();
- writeImage<MaskImageType>(skeletonFilter->GetOutput(), "skel.mhd");
- writeImage<MaskImageType>(skeletonFilter->GetThinning(), "skel2.mhd");
- */
-
- //-----------------------------------------------------
- StartNewStep("Left limits with bronchus (slice by slice)");
- // Select LeftLabel (set right label to 0)
- MaskImagePointer temp = SetBackground<MaskImageType, MaskImageType>(m_working_trachea, m_working_trachea, rightLabel, 0);
- writeImage<MaskImageType>(temp, "temp1.mhd");
-
- m_working_image =
- clitk::SliceBySliceRelativePosition<MaskImageType>(m_working_image,
- temp,
- 2, GetFuzzyThreshold(), "RightTo");
- /*
- typedef clitk::SliceBySliceRelativePositionFilter<MaskImageType> SliceRelPosFilterType;
- typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New();
- sliceRelPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- sliceRelPosFilter->VerboseStepOn();
- sliceRelPosFilter->WriteStepOn();
- sliceRelPosFilter->SetInput(m_working_image);
- sliceRelPosFilter->SetInputObject(temp);
- sliceRelPosFilter->SetDirection(2);
- sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThreshold());
- sliceRelPosFilter->SetOrientationTypeString("RightTo");
- sliceRelPosFilter->Update();
- m_working_image = sliceRelPosFilter->GetOutput();*/
- writeImage<MaskImageType>(m_working_image, "afterslicebyslice.mhd");
-
-
- typedef clitk::SliceBySliceRelativePositionFilter<MaskImageType> SliceRelPosFilterType;
- typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New();
-
-
- //-----------------------------------------------------
- StartNewStep("Right limits with bronchus (slice by slice)");
- // Select LeftLabel (set right label to 0)
- temp = SetBackground<MaskImageType, MaskImageType>(m_working_trachea, m_working_trachea, leftLabel, 0);
- writeImage<MaskImageType>(temp, "temp2.mhd");
-
- sliceRelPosFilter = SliceRelPosFilterType::New();
- sliceRelPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- sliceRelPosFilter->VerboseStepOff();
- sliceRelPosFilter->WriteStepOff();
- sliceRelPosFilter->SetInput(m_working_image);
- sliceRelPosFilter->SetInputObject(temp);
- sliceRelPosFilter->SetDirection(2);
- sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThreshold());
- sliceRelPosFilter->SetOrientationTypeString("LeftTo");
- sliceRelPosFilter->Update();
- m_working_image = sliceRelPosFilter->GetOutput();
- writeImage<MaskImageType>(m_working_image, "afterslicebyslice.mhd");
- DD("end");
- m_output = m_working_image;
- StopCurrentStep<MaskImageType>(m_output);
-
- //-----------------------------------------------------
- StartNewStep("Trial Post position according to trachea");
- // Post: do not extend past the post wall of the 2 main bronchi
- sliceRelPosFilter = SliceRelPosFilterType::New();
- sliceRelPosFilter->SetCurrentStepBaseId(this->GetCurrentStepId());
- sliceRelPosFilter->VerboseStepOn();
- sliceRelPosFilter->WriteStepOn();
- sliceRelPosFilter->SetInput(m_output);
- sliceRelPosFilter->SetInputObject(m_trachea);
- sliceRelPosFilter->SetDirection(2);
- sliceRelPosFilter->SetFuzzyThreshold(GetFuzzyThreshold());
- // It means "not PostTo" (different from AntTo)
- sliceRelPosFilter->NotFlagOn();
- //sliceRelPosFilter->SetOrientationType(SliceRelPosFilterType::RelPosFilterType::PostTo);
- sliceRelPosFilter->SetOrientationTypeString("PostTo");
- sliceRelPosFilter->Update();
- m_output = sliceRelPosFilter->GetOutput();
- writeImage<MaskImageType>(m_output, "postrelpos.mhd");
-
-
- // Set output image information (required)
- MaskImagePointer outputImage = this->GetOutput(0);
- outputImage->SetRegions(m_working_image->GetLargestPossibleRegion());
- outputImage->SetOrigin(m_working_image->GetOrigin());
- outputImage->SetRequestedRegion(m_working_image->GetLargestPossibleRegion());
- DD("end2");
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class TImageType>
-void
-clitk::ExtractLymphStationsFilter<TImageType>::
-GenerateInputRequestedRegion() {
- DD("GenerateInputRequestedRegion");
- // Call default
- // Superclass::GenerateInputRequestedRegion();
- // Following needed because output region can be greater than input (trachea)
- //ImagePointer mediastinum = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(0));
- //ImagePointer trachea = dynamic_cast<TImageType*>(itk::ProcessObject::GetInput(1));
- DD("set reg");
- m_mediastinum->SetRequestedRegion(m_mediastinum->GetLargestPossibleRegion());
- m_trachea->SetRequestedRegion(m_trachea->GetLargestPossibleRegion());
- DD("end");
-}
-//--------------------------------------------------------------------
-
-
-//--------------------------------------------------------------------
-template <class TImageType>
-void
-clitk::ExtractLymphStationsFilter<TImageType>::
-GenerateData() {
- DD("GenerateData");
- // Final Step -> graft output (if SetNthOutput => redo)
- this->GraftOutput(m_output);
-}
-//--------------------------------------------------------------------
-
-
-#endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX
f->SetWriteStepFlag(mArgsInfo.writeStep_flag);
f->SetVerboseMemoryFlag(mArgsInfo.verboseMemory_flag);
f->SetAFDBFilename(mArgsInfo.afdb_arg);
-
- f->SetComputeStationsSupportsFlag(!mArgsInfo.nosupport_flag);
+ if (mArgsInfo.afdb_path_given) f->SetAFDBPath(mArgsInfo.afdb_path_arg);
+ f->SetForceSupportsFlag(mArgsInfo.force_support_flag);
+ f->SetSupportLimitsFilename(mArgsInfo.support_limits_arg);
+ f->SetCheckSupportFlag(mArgsInfo.check_support_limits_flag);
// Station 8
//f->SetDistanceMaxToAnteriorPartOfTheSpine(mArgsInfo.S8_maxAntSpine_arg);
// Generic RelativePosition processes
output = this->ApplyRelativePositionList("Mediastinum", output);
+
//--------------------------------------------------------------------
+ // FIXME --> do not put this limits here !
+ /*
// Step : SI limits It is better to do this limit *AFTER* the
// RelativePosition to avoid some issue due to superior boundaries.
this->StartNewStep("[Mediastinum] Keep inferior to CricoidCartilag");
}
output = clitk::CropImageRemoveGreaterThan<MaskImageType>(output, 2, p[2], true, this->GetBackgroundValue());
this->template StopCurrentStep<MaskImageType>(output);
+ */
//--------------------------------------------------------------------
// Step: Get CCL
filter->SetArgsInfo(args_info);
- try {
- filter->Update();
- } catch(std::runtime_error e) {
- std::cout << e.what() << std::endl;
- }
+ CLITK_TRY_CATCH_EXIT(filter->Update());
return EXIT_SUCCESS;
} // This is the end, my friend
//-----------------------------------------------------------
// Constructor
//-----------------------------------------------------------
- ExtractPatientGenericFilter::ExtractPatientGenericFilter()
+ ExtractPatientGenericFilter:::ExtractPatientGenericFilter()
{
m_Verbose=false;
m_InputFileName="";
{
m_AFDB = NULL;
SetAFDBFilename("default.afdb");
+ SetAFDBPath("./");
}
//--------------------------------------------------------------------
void clitk::FilterWithAnatomicalFeatureDatabaseManagement::LoadAFDB()
{
GetAFDB()->SetFilename(GetAFDBFilename());
+ GetAFDB()->SetPath(GetAFDBPath());
try {
GetAFDB()->Load();
} catch (clitk::ExceptionObject e) {
//--------------------------------------------------------------------
-clitk::AnatomicalFeatureDatabase * clitk::FilterWithAnatomicalFeatureDatabaseManagement::GetAFDB()
+clitk::AnatomicalFeatureDatabase::Pointer clitk::FilterWithAnatomicalFeatureDatabaseManagement::GetAFDB()
{
- if (m_AFDB == NULL) {
- m_AFDB = new clitk::AnatomicalFeatureDatabase;
+ if (!m_AFDB) {
+ m_AFDB = clitk::AnatomicalFeatureDatabase::New();
}
return m_AFDB;
}
itkTypeMacro(FilterWithAnatomicalFeatureDatabaseManagement, Object);
// Set/Get filename
- itkBooleanMacro(AFDBFilenameGivenFlag);
- itkSetMacro(AFDBFilenameGivenFlag, bool);
- itkGetConstMacro(AFDBFilenameGivenFlag, bool);
- GGO_DefineOption_Flag(afdb, SetAFDBFilenameGivenFlag);
+ // itkBooleanMacro(AFDBFilenameGivenFlag);
+ // itkSetMacro(AFDBFilenameGivenFlag, bool);
+ // itkGetConstMacro(AFDBFilenameGivenFlag, bool);
+ // GGO_DefineOption_Flag(afdb, SetAFDBFilenameGivenFlag);
+
+ // itkBooleanMacro(AFDBPathGivenFlag);
+ // itkSetMacro(AFDBPathGivenFlag, bool);
+ // itkGetConstMacro(AFDBPathGivenFlag, bool);
+ // GGO_DefineOption_Flag(afdb_path, SetAFDBPathGivenFlag);
itkSetMacro(AFDBFilename, std::string);
itkGetConstMacro(AFDBFilename, std::string);
- GGO_DefineOption_WithTest(afdb, SetAFDBFilename, std::string, AFDBFilenameGivenFlag);
+ // GGO_DefineOption_WithTest(afdb, SetAFDBFilename, std::string, AFDBFilenameGivenFlag);
+
+ itkSetMacro(AFDBPath, std::string);
+ itkGetConstMacro(AFDBPath, std::string);
+ // GGO_DefineOption_WithTest(afdb_path, SetAFDBPath, std::string, AFDBPathGivenFlag);
void WriteAFDB();
void LoadAFDB();
- AnatomicalFeatureDatabase * GetAFDB();
+ AnatomicalFeatureDatabase::Pointer GetAFDB();
void SetAFDB(AnatomicalFeatureDatabase * a) { m_AFDB = a; }
protected:
virtual ~FilterWithAnatomicalFeatureDatabaseManagement() {}
std::string m_AFDBFilename;
- bool m_AFDBFilenameGivenFlag;
- clitk::AnatomicalFeatureDatabase * m_AFDB;
+ // bool m_AFDBFilenameGivenFlag;
+ std::string m_AFDBPath;
+ // bool m_AFDBPathGivenFlag;
+ clitk::AnatomicalFeatureDatabase::Pointer m_AFDB;
private:
FilterWithAnatomicalFeatureDatabaseManagement(const Self&); //purposely not implemented
// Passing the arguments to the generic filter
filter->SetArgsInfo(args_info);
- try {
- filter->Update();
- } catch(std::runtime_error e) {
- std::cout << e.what() << std::endl;
- }
+ CLITK_TRY_CATCH_EXIT(filter->Update());
return EXIT_SUCCESS;
} // This is the end, my friend
clitk::MotionMaskGenericFilter::Pointer genericFilter=clitk::MotionMaskGenericFilter::New();
genericFilter->SetArgsInfo(args_info);
- genericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(genericFilter->Update());
return EXIT_SUCCESS;
}// end main
// Loop on RelativePositionList of operations
std::string s = GetInputName();
for(uint i=0; i<mArgsInfoList.size(); i++) {
- std::string text = "["+s+"] limits "+
- mArgsInfoList[i].orientation_arg[0]+" "+
- mArgsInfoList[i].object_arg+" "+
- toString(mArgsInfoList[i].threshold_arg);
+ std::string text = "["+s+"] limits ";
+ if (mArgsInfoList[i].orientation_given) text += std::string(mArgsInfoList[i].orientation_arg[0])+" ";
+ else text = text+"("+toString(mArgsInfoList[i].angle1_arg)+" "+
+ toString(mArgsInfoList[i].angle2_arg)+" "+
+ (mArgsInfoList[i].inverse_flag==true?"true":"false")+") ";
+ text = text+mArgsInfoList[i].object_arg+" "+toString(mArgsInfoList[i].threshold_arg);
if (mArgsInfoList[i].sliceBySlice_flag) {
text += " slice by slice";
}
clitk::RelativePositionList<TImageType>::
SetFilterOptions(typename RelPosFilterType::Pointer filter, ArgsInfoType & options) {
- if (options.orientation_given != 1) {
- DD("ERRROR DEBUG TODO no more than 1 orientation yet");
- exit(0);
+ if (options.orientation_given > 1) {
+ clitkExceptionMacro("Error in the RelPos options. I need a single --orientation (for the moment)."
+ << " Current options are for obj = '" << options.object_arg
+ << "', threshold = " << options.threshold_arg << std::endl);
}
ImagePointer object = GetAFDB()->template GetImage<ImageType>(options.object_arg);
filter->SetVerboseImageSizeFlag(GetVerboseImageSizeFlag());
filter->SetFuzzyThreshold(options.threshold_arg);
filter->SetInverseOrientationFlag(options.inverse_flag); // MUST BE BEFORE AddOrientationTypeString
- for(uint i=0; i<options.orientation_given; i++)
- filter->AddOrientationTypeString(options.orientation_arg[i]);
+
+ if (options.orientation_given == 1) {
+ for(uint i=0; i<options.orientation_given; i++)
+ filter->AddOrientationTypeString(options.orientation_arg[i]);
+ }
+ else {
+ if (options.angle1_given && options.angle2_given) {
+ filter->AddAnglesInDeg(options.angle1_arg, options.angle2_arg);
+ }
+ else {
+ clitkExceptionMacro("Error in the RelPos options. I need --orientation or (--angle1 and --angle2)."
+ << " Current options are for obj = '" << options.object_arg
+ << "', threshold = " << options.threshold_arg << std::endl);
+ }
+ }
filter->SetIntermediateSpacing(options.spacing_arg);
if (options.spacing_arg == -1) filter->IntermediateSpacingFlagOff();
else filter->IntermediateSpacingFlagOn();
TARGET_LINK_LIBRARIES(clitkRelativePosition clitkCommon ${ITK_LIBRARIES})
SET(TOOLS_INSTALL ${TOOLS_INSTALL} clitkRelativePosition)
+ WRAP_GGO(clitkLabelImageOverlapMeasure_GGO_C clitkLabelImageOverlapMeasure.ggo)
+ ADD_EXECUTABLE(clitkLabelImageOverlapMeasure clitkLabelImageOverlapMeasure.cxx ${clitkLabelImageOverlapMeasure_GGO_C})
+ TARGET_LINK_LIBRARIES(clitkLabelImageOverlapMeasure clitkSegmentationGgoLib clitkCommon ${ITK_LIBRARIES} )
+ SET(TOOLS_INSTALL ${TOOLS_INSTALL} clitkLabelImageOverlapMeasure)
+
+ WRAP_GGO(clitkRelativePositionAnalyzer_GGO_C clitkRelativePositionAnalyzer.ggo)
+ ADD_EXECUTABLE(clitkRelativePositionAnalyzer clitkRelativePositionAnalyzer.cxx ${clitkRelativePositionAnalyzer_GGO_C})
+ TARGET_LINK_LIBRARIES(clitkRelativePositionAnalyzer clitkSegmentationGgoLib clitkCommon ${ITK_LIBRARIES})
+ SET(TOOLS_INSTALL ${TOOLS_INSTALL} clitkRelativePositionAnalyzer)
+
+ WRAP_GGO(clitkRelativePositionDataBaseAnalyzer_GGO_C clitkRelativePositionDataBaseAnalyzer.ggo)
+ ADD_EXECUTABLE(clitkRelativePositionDataBaseAnalyzer ../itk/clitkRelativePositionDataBase.cxx clitkRelativePositionDataBaseAnalyzer.cxx ${clitkRelativePositionDataBaseAnalyzer_GGO_C})
+ TARGET_LINK_LIBRARIES(clitkRelativePositionDataBaseAnalyzer clitkSegmentationGgoLib clitkCommon ${ITK_LIBRARIES})
+ SET(TOOLS_INSTALL ${TOOLS_INSTALL} clitkRelativePositionDataBaseAnalyzer)
+
+ WRAP_GGO(clitkRelativePositionDataBaseBuilder_GGO_C clitkRelativePositionDataBaseBuilder.ggo)
+ ADD_EXECUTABLE(clitkRelativePositionDataBaseBuilder clitkRelativePositionDataBaseBuilder.cxx ${clitkRelativePositionDataBaseBuilder_GGO_C})
+ TARGET_LINK_LIBRARIES(clitkRelativePositionDataBaseBuilder clitkSegmentationGgoLib clitkCommon ${ITK_LIBRARIES})
+ SET(TOOLS_INSTALL ${TOOLS_INSTALL} clitkRelativePositionDataBaseBuilder)
+
WRAP_GGO(clitkTransformLandmarks_GGO_C clitkTransformLandmarks.ggo)
ADD_EXECUTABLE(clitkTransformLandmarks clitkTransformLandmarks.cxx ${clitkTransformLandmarks_GGO_C})
TARGET_LINK_LIBRARIES(clitkTransformLandmarks clitkCommon ${ITK_LIBRARIES})
FilterType::Pointer genericFilter = FilterType::New();
genericFilter->SetArgsInfo(args_info);
- genericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(genericFilter->Update());
return EXIT_SUCCESS;
}// end main
filter->SetArgsInfo(args_info);
- try {
- filter->Update();
- } catch(std::runtime_error e) {
- std::cout << e.what() << std::endl;
- }
+ CLITK_TRY_CATCH_EXIT(filter->Update());
return EXIT_SUCCESS;
} // This is the end, my friend
clitk::CombineImageGenericFilter::Pointer genericFilter=clitk::CombineImageGenericFilter::New();
genericFilter->SetArgsInfo(args_info);
- genericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(genericFilter->Update());
return EXIT_SUCCESS;
}// end main
//JV how to pass for different dims?
//ComposeVFGenericFilter->SetEdgePaddingValue(args_info.pad_arg);
ComposeVFGenericFilter->SetVerbose(args_info.verbose_flag);
- ComposeVFGenericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(ComposeVFGenericFilter->Update());
return EXIT_SUCCESS;
}
- BSD See included LICENSE.txt file
- CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
-===========================================================================**/
+ ===========================================================================**/
#ifndef clitkGetSpacingGenericFilter_cxx
#define clitkGetSpacingGenericFilter_cxx
-/* =================================================
- * @file clitkGetSpacingGenericFilter.cxx
- * @author
- * @date
- *
- * @brief
- *
- ===================================================*/
-
#include "clitkGetSpacingGenericFilter.h"
-
namespace clitk
{
-
//-----------------------------------------------------------
// Constructor
//-----------------------------------------------------------
// Go !
filter->SetArgsInfo(args_info);
- filter->Update();
+ CLITK_TRY_CATCH_EXIT(filter->Update());
// this is the end my friend
return EXIT_SUCCESS;
if (args_info.type_given) filter->SetOutputPixelType(args_info.type_arg);
// Go !
- filter->Update();
+ CLITK_TRY_CATCH_EXIT(filter->Update());
// this is the end my friend
return 0;
FilterType::Pointer genericFilter = FilterType::New();
genericFilter->SetArgsInfo(args_info);
- genericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(genericFilter->Update());
return EXIT_SUCCESS;
}// end main
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+
+// clitk
+#include "clitkIO.h"
+#include "clitkLabelImageOverlapMeasure_ggo.h"
+#include "clitkLabelImageOverlapMeasureGenericFilter.h"
+
+//--------------------------------------------------------------------
+int main(int argc, char * argv[]) {
+
+ // Init command line
+ GGO(clitkLabelImageOverlapMeasure, args_info);
+ CLITK_INIT;
+
+ // Filter
+ typedef clitk::LabelImageOverlapMeasureGenericFilter<args_info_clitkLabelImageOverlapMeasure> FilterType;
+ FilterType::Pointer filter = FilterType::New();
+
+ filter->SetArgsInfo(args_info);
+
+ try {
+ filter->Update();
+ } catch(std::runtime_error e) {
+ std::cout << e.what() << std::endl;
+ return EXIT_FAILURE;
+ }
+
+ return EXIT_SUCCESS;
+} // This is the end, my friend
+//--------------------------------------------------------------------
--- /dev/null
+#File clitkLabelImageOverlapMeasure.ggo
+package "clitkLabelImageOverlapMeasure"
+version "1.0"
+purpose "Compute Dice and other coefficients between label images"
+
+section "General options"
+option "config" - "Config file" string no
+option "verbose" v "Verbose" flag off
+option "imagetypes" - "Display allowed image types" flag off
+
+section "Input"
+option "input1" i "Input mask 1" string yes
+option "input2" j "Input mask 2" string yes
+
+option "label1" l "Label in input1" int no default="1"
+option "label2" m "Label in input2" int no default="1"
+option "BG" b "Background value" int no default="0"
+
+
+
+
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+
+#ifndef CLITKLABELIMAGEOVERLAPMEASUREGENERICFILTER_H
+#define CLITKLABELIMAGEOVERLAPMEASUREGENERICFILTER_H
+
+// clitk
+#include "clitkImageToImageGenericFilter.h"
+#include "clitkLabelImageOverlapMeasureFilter.h"
+#include "clitkBoundingBoxUtils.h"
+#include "clitkCropLikeImageFilter.h"
+
+//--------------------------------------------------------------------
+namespace clitk
+{
+
+ template<class ArgsInfoType>
+ class ITK_EXPORT LabelImageOverlapMeasureGenericFilter:
+ public ImageToImageGenericFilter<LabelImageOverlapMeasureGenericFilter<ArgsInfoType> >
+ {
+ public:
+ //--------------------------------------------------------------------
+ LabelImageOverlapMeasureGenericFilter();
+
+ //--------------------------------------------------------------------
+ typedef ImageToImageGenericFilter<LabelImageOverlapMeasureGenericFilter<ArgsInfoType> > Superclass;
+ typedef LabelImageOverlapMeasureGenericFilter Self;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ //--------------------------------------------------------------------
+ itkNewMacro(Self);
+ itkTypeMacro(LabelImageOverlapMeasureGenericFilter, LightObject);
+
+ //--------------------------------------------------------------------
+ void SetArgsInfo(const ArgsInfoType & a);
+ template<class FilterType>
+ void SetOptionsFromArgsInfoToFilter(FilterType * f) ;
+
+ //--------------------------------------------------------------------
+ // Main function called each time the filter is updated
+ template<class ImageType>
+ void UpdateWithInputImageType();
+
+ protected:
+ template<unsigned int Dim> void InitializeImageType();
+ ArgsInfoType mArgsInfo;
+
+ private:
+ LabelImageOverlapMeasureGenericFilter(const Self&); //purposely not implemented
+ void operator=(const Self&); //purposely not implemented
+
+ };// end class
+ //--------------------------------------------------------------------
+} // end namespace clitk
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include "clitkLabelImageOverlapMeasureGenericFilter.txx"
+#endif
+
+#endif // #define CLITKRELATIVEPOSITIONANALYZERGENERICFILTER_H
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+clitk::LabelImageOverlapMeasureGenericFilter<ArgsInfoType>::
+LabelImageOverlapMeasureGenericFilter():
+ ImageToImageGenericFilter<Self>("LabelImageOverlapMeasure")
+{
+ // Default values
+ cmdline_parser_clitkLabelImageOverlapMeasure_init(&mArgsInfo);
+ //InitializeImageType<2>();
+ InitializeImageType<3>();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<unsigned int Dim>
+void clitk::LabelImageOverlapMeasureGenericFilter<ArgsInfoType>::
+InitializeImageType()
+{
+ ADD_IMAGE_TYPE(Dim, uchar);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+void clitk::LabelImageOverlapMeasureGenericFilter<ArgsInfoType>::
+SetArgsInfo(const ArgsInfoType & a)
+{
+ mArgsInfo=a;
+ SetIOVerbose(mArgsInfo.verbose_flag);
+ if (mArgsInfo.imagetypes_flag) this->PrintAvailableImageTypes();
+ if (mArgsInfo.input1_given) AddInputFilename(mArgsInfo.input1_arg);
+ if (mArgsInfo.input2_given) AddInputFilename(mArgsInfo.input2_arg);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+// Update with the number of dimensions and the pixeltype
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class FilterType>
+void clitk::LabelImageOverlapMeasureGenericFilter<ArgsInfoType>::
+SetOptionsFromArgsInfoToFilter(FilterType * f)
+{
+ f->SetLabel1(mArgsInfo.label1_arg);
+ f->SetLabel2(mArgsInfo.label2_arg);
+}
+
+//--------------------------------------------------------------------
+// Update with the number of dimensions and the pixeltype
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class ImageType>
+void clitk::LabelImageOverlapMeasureGenericFilter<ArgsInfoType>::
+UpdateWithInputImageType()
+{
+ // Reading input
+ typename ImageType::Pointer input1 = this->template GetInput<ImageType>(0);
+ typename ImageType::Pointer input2 = this->template GetInput<ImageType>(1);
+
+ // Create filter
+ typedef clitk::LabelImageOverlapMeasureFilter<ImageType> FilterType;
+ typename FilterType::Pointer filter = FilterType::New();
+
+ // Set global Options
+ filter->SetInput(0, input1);
+ filter->SetInput(1, input2);
+ SetOptionsFromArgsInfoToFilter<FilterType>(filter);
+
+ // Go !
+ filter->Update();
+
+ // Write/Save results
+ // typename ImageType::Pointer output = filter->GetOutput();
+ // this->template SetNextOutput<ImageType>(output);
+}
+//--------------------------------------------------------------------
+
+
FilterType::Pointer genericFilter = FilterType::New();
genericFilter->SetArgsInfo(args_info);
- genericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(genericFilter->Update());
return EXIT_SUCCESS;
}// end main
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+
+// clitk
+#include "clitkRelativePositionAnalyzer_ggo.h"
+#include "clitkRelativePositionAnalyzerGenericFilter.h"
+
+//--------------------------------------------------------------------
+int main(int argc, char * argv[]) {
+
+ // Init command line
+ GGO(clitkRelativePositionAnalyzer, args_info);
+ CLITK_INIT;
+
+ // Filter
+ typedef clitk::RelativePositionAnalyzerGenericFilter<args_info_clitkRelativePositionAnalyzer> FilterType;
+ FilterType::Pointer filter = FilterType::New();
+ filter->SetArgsInfo(args_info);
+
+ // Go !
+ try {
+ filter->Update();
+ } catch(std::runtime_error e) {
+ std::cout << e.what() << std::endl;
+ return EXIT_FAILURE;
+ }
+
+ return EXIT_SUCCESS;
+} // This is the end, my friend
+//--------------------------------------------------------------------
--- /dev/null
+#File clitkRelativePositionAnalyzer.ggo
+package "clitkRelativePositionAnalyzer"
+version "1.0"
+purpose "Analyze relative position of a target according to structures"
+
+section "General options"
+option "config" - "Config file" string no
+option "verbose" v "Verbose" flag off
+option "imagetypes" - "Display allowed image types" flag off
+
+section "Input/Output"
+option "support" i "Input mask support" string yes
+option "object" j "Input mask object" string yes
+option "target" t "Input mask target" string yes
+option "output" o "Output image " string yes
+
+section "Options for building the relative positions"
+option "bins" - "Number of histo bins for fuzzy map" int default="100" no
+option "tol" - "Target area loss tolerance (|0-1])" double default="0.01" no
+
+option "angle1" a "Angle 1 (deg)" double no default="0"
+option "angle2" b "Angle 2 (deg)" double no default="0"
+option "inverse" n "Not flag : inverse of the orientation" flag off
+
+
+section "Compute some statistic on afdb"
+option "afdb" - "Input Anatomical Feature DB" string no multiple
+option "afdb_path" - "Path to search image in afdb" string no multiple
+option "station" - "Station name" string no
+
+
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+
+#ifndef CLITKRELATIVEPOSITIONANALYZERGENERICFILTER_H
+#define CLITKRELATIVEPOSITIONANALYZERGENERICFILTER_H
+
+// clitk
+#include "clitkIO.h"
+#include "clitkImageToImageGenericFilter.h"
+#include "clitkRelativePositionAnalyzerFilter.h"
+#include "clitkSliceBySliceRelativePositionFilter.h"
+
+//--------------------------------------------------------------------
+namespace clitk
+{
+
+ template<class ArgsInfoType>
+ class ITK_EXPORT RelativePositionAnalyzerGenericFilter:
+ public ImageToImageGenericFilter<RelativePositionAnalyzerGenericFilter<ArgsInfoType> >
+ {
+ public:
+ //--------------------------------------------------------------------
+ RelativePositionAnalyzerGenericFilter();
+ ~RelativePositionAnalyzerGenericFilter() {}
+
+ //--------------------------------------------------------------------
+ typedef ImageToImageGenericFilter<RelativePositionAnalyzerGenericFilter<ArgsInfoType> > Superclass;
+ typedef RelativePositionAnalyzerGenericFilter Self;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ //--------------------------------------------------------------------
+ itkNewMacro(Self);
+ itkTypeMacro(RelativePositionAnalyzerGenericFilter, LightObject);
+
+ //--------------------------------------------------------------------
+ void SetArgsInfo(const ArgsInfoType & a);
+ template<class FilterType>
+ void SetOptionsFromArgsInfoToFilter(FilterType * f) ;
+
+ //--------------------------------------------------------------------
+ // Main function called each time the filter is updated
+ template<class ImageType>
+ void UpdateWithInputImageType();
+
+ protected:
+ template<unsigned int Dim> void InitializeImageType();
+ ArgsInfoType mArgsInfo;
+
+ private:
+ RelativePositionAnalyzerGenericFilter(const Self&); //purposely not implemented
+ void operator=(const Self&); //purposely not implemented
+
+ };// end class
+ //--------------------------------------------------------------------
+} // end namespace clitk
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include "clitkRelativePositionAnalyzerGenericFilter.txx"
+#endif
+
+#endif // #define CLITKRELATIVEPOSITIONANALYZERGENERICFILTER_H
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+clitk::RelativePositionAnalyzerGenericFilter<ArgsInfoType>::
+RelativePositionAnalyzerGenericFilter():
+ ImageToImageGenericFilter<Self>("RelativePositionAnalyzer")
+{
+ // Default values
+ cmdline_parser_clitkRelativePositionAnalyzer_init(&mArgsInfo);
+ InitializeImageType<3>();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<unsigned int Dim>
+void clitk::RelativePositionAnalyzerGenericFilter<ArgsInfoType>::
+InitializeImageType()
+{
+ ADD_IMAGE_TYPE(Dim, uchar);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+void clitk::RelativePositionAnalyzerGenericFilter<ArgsInfoType>::
+SetArgsInfo(const ArgsInfoType & a)
+{
+ mArgsInfo=a;
+ SetIOVerbose(mArgsInfo.verbose_flag);
+ if (mArgsInfo.imagetypes_flag) this->PrintAvailableImageTypes();
+ if (mArgsInfo.support_given) AddInputFilename(mArgsInfo.support_arg);
+ if (mArgsInfo.object_given) AddInputFilename(mArgsInfo.object_arg);
+ if (mArgsInfo.target_given) AddInputFilename(mArgsInfo.target_arg);
+ if (mArgsInfo.output_given) AddOutputFilename(mArgsInfo.output_arg);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+// Update with the number of dimensions and the pixeltype
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class FilterType>
+void clitk::RelativePositionAnalyzerGenericFilter<ArgsInfoType>::
+SetOptionsFromArgsInfoToFilter(FilterType * f)
+{
+ f->SetNumberOfBins(mArgsInfo.bins_arg);
+ f->SetAreaLossTolerance(mArgsInfo.tol_arg);
+ clitk::RelativePositionDirectionType d;
+ if (mArgsInfo.angle1_given) d.angle1 = clitk::deg2rad(mArgsInfo.angle1_arg);
+ if (mArgsInfo.angle2_given) d.angle2 = clitk::deg2rad(mArgsInfo.angle2_arg);
+ if (mArgsInfo.inverse_given) d.notFlag = clitk::deg2rad(mArgsInfo.inverse_flag);
+ f->SetDirection(d);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+// Update with the number of dimensions and the pixeltype
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class ImageType>
+void clitk::RelativePositionAnalyzerGenericFilter<ArgsInfoType>::
+UpdateWithInputImageType()
+{
+ // Create filter
+ typedef clitk::RelativePositionAnalyzerFilter<ImageType> FilterType;
+ typename FilterType::Pointer filter = FilterType::New();
+
+ // Reading input
+ typename ImageType::Pointer support = this->template GetInput<ImageType>(0);
+ typename ImageType::Pointer object = this->template GetInput<ImageType>(1);
+ typename ImageType::Pointer target = this->template GetInput<ImageType>(2);
+ filter->SetInputSupport(support);
+ filter->SetInputObject(object);
+ filter->SetInputTarget(target);
+
+ // Set global Options
+ SetOptionsFromArgsInfoToFilter<FilterType>(filter);
+
+ // Go !
+ filter->Update();
+
+ // Display output
+ filter->GetInfo().Println();
+ filter->GetInfoReverse().Println();
+}
+//--------------------------------------------------------------------
+
+
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+
+// clitk
+#include "clitkRelativePositionDataBaseAnalyzer_ggo.h"
+#include "clitkIO.h"
+#include "clitkRelativePositionDataBaseAnalyzerFilter.h"
+
+//--------------------------------------------------------------------
+int main(int argc, char * argv[]) {
+
+ // Init command line
+ GGO(clitkRelativePositionDataBaseAnalyzer, args_info);
+ CLITK_INIT;
+
+ // Filter
+ typedef clitk::RelativePositionDataBaseAnalyzerFilter FilterType;
+ FilterType filter;
+ filter.SetDatabaseFilename(args_info.db_arg);
+ filter.SetStationName(args_info.station_arg);
+ // filter.SetStationName(args_info.object_arg);//FIXME
+ filter.SetOutputFilename(args_info.output_arg);
+
+ try {
+ filter.Update();
+ } catch(std::runtime_error e) {
+ std::cout << e.what() << std::endl;
+ return EXIT_FAILURE;
+ }
+
+ return EXIT_SUCCESS;
+} // This is the end, my friend
+//--------------------------------------------------------------------
--- /dev/null
+#File clitkRelativePositionDataBaseAnalyzer.ggo
+package "clitkRelativePositionDataBaseAnalyzer"
+version "1.0"
+purpose "Analyze a DB of relative positions for several patient and structures"
+
+section "General options"
+option "config" - "Config file" string no
+option "verbose" v "Verbose" flag off
+
+section "Input/Output"
+option "db" i "Input database of RelPost" string yes
+option "station" s "Studied station name" string yes
+option "object" - "Studied object/struct name" string no
+option "output" o "Output rpl file" string yes
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+
+// clitk
+#include "clitkRelativePositionDataBaseBuilder_ggo.h"
+#include "clitkRelativePositionDataBaseBuilderGenericFilter.h"
+
+//--------------------------------------------------------------------
+int main(int argc, char * argv[]) {
+
+ // Init command line
+ GGO(clitkRelativePositionDataBaseBuilder, args_info);
+ CLITK_INIT;
+
+ // Filter
+ typedef clitk::RelativePositionDataBaseBuilderGenericFilter<args_info_clitkRelativePositionDataBaseBuilder> FilterType;
+ FilterType::Pointer filter = FilterType::New();
+
+ // Set options
+ filter->SetArgsInfo(args_info);
+
+ // Add an input to determine the type of image
+ NewAFDB(afdb, args_info.afdb_arg);
+ std::string f = afdb->GetTagValue(args_info.supportName_arg);
+ f = std::string(args_info.afdb_path_arg)+"/"+f;
+ filter->AddInputFilename(f);
+
+ try {
+ filter->Update();
+ } catch(std::runtime_error e) {
+ std::cout << e.what() << std::endl;
+ return EXIT_FAILURE;
+ }
+
+ return EXIT_SUCCESS;
+} // This is the end, my friend
+//--------------------------------------------------------------------
--- /dev/null
+#File clitkRelativePositionDataBaseBuilder.ggo
+package "clitkRelativePositionDataBaseBuilder"
+version "1.0"
+purpose "Analyze relative position of a target according to structures"
+
+section "General options"
+option "config" - "Config file" string no
+option "verbose" v "Verbose" flag off
+option "imagetypes" - "Display allowed image types" flag off
+
+section "Input/Output"
+option "afdb" a "Input Anatomical Feature DB" string yes
+option "afdb_path" - "Path to search image in afdb" string no
+
+option "supportName" i "Input mask support name in afdb" string yes
+option "targetName" t "Input mask target name in afdb" string yes
+option "objectName" j "Input mask object name in afdb" string yes multiple
+
+option "output" o "Output image " string yes
+
+section "Options for building the relative positions"
+option "bins" b "Number of histo bins for fuzzy map" int default="100" no
+option "nb" n "Number of angles to test" int default="4" no
+option "tol" - "Target area loss tolerance (|0-1])" double default="0.01" no
+
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+
+#ifndef CLITKRelativePositionDataBaseBuilderGENERICFILTER_H
+#define CLITKRelativePositionDataBaseBuilderGENERICFILTER_H
+
+// clitk
+#include "clitkIO.h"
+#include "clitkImageToImageGenericFilter.h"
+#include "clitkRelativePositionDataBaseBuilderFilter.h"
+#include "clitkSliceBySliceRelativePositionFilter.h"
+
+//--------------------------------------------------------------------
+namespace clitk
+{
+
+ template<class ArgsInfoType>
+ class ITK_EXPORT RelativePositionDataBaseBuilderGenericFilter:
+ public ImageToImageGenericFilter<RelativePositionDataBaseBuilderGenericFilter<ArgsInfoType> >
+ {
+ public:
+ //--------------------------------------------------------------------
+ RelativePositionDataBaseBuilderGenericFilter();
+
+ //--------------------------------------------------------------------
+ typedef ImageToImageGenericFilter<RelativePositionDataBaseBuilderGenericFilter<ArgsInfoType> > Superclass;
+ typedef RelativePositionDataBaseBuilderGenericFilter Self;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ //--------------------------------------------------------------------
+ itkNewMacro(Self);
+ itkTypeMacro(RelativePositionDataBaseBuilderGenericFilter, LightObject);
+
+ //--------------------------------------------------------------------
+ void SetArgsInfo(const ArgsInfoType & a);
+ template<class FilterType>
+ void SetOptionsFromArgsInfoToFilter(FilterType * f) ;
+
+ //--------------------------------------------------------------------
+ // Main function called each time the filter is updated
+ template<class ImageType>
+ void UpdateWithInputImageType();
+
+ protected:
+ template<unsigned int Dim> void InitializeImageType();
+ ArgsInfoType mArgsInfo;
+
+ private:
+ RelativePositionDataBaseBuilderGenericFilter(const Self&); //purposely not implemented
+ void operator=(const Self&); //purposely not implemented
+
+ };// end class
+ //--------------------------------------------------------------------
+} // end namespace clitk
+
+#ifndef ITK_MANUAL_INSTANTIATION
+#include "clitkRelativePositionDataBaseBuilderGenericFilter.txx"
+#endif
+
+#endif // #define CLITKRelativePositionDataBaseBuilderGENERICFILTER_H
--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+ ===========================================================================**/
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+clitk::RelativePositionDataBaseBuilderGenericFilter<ArgsInfoType>::
+RelativePositionDataBaseBuilderGenericFilter():
+ ImageToImageGenericFilter<Self>("RelativePositionDataBaseBuilder")
+{
+ // Default values
+ cmdline_parser_clitkRelativePositionDataBaseBuilder_init(&mArgsInfo);
+ InitializeImageType<3>();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<unsigned int Dim>
+void clitk::RelativePositionDataBaseBuilderGenericFilter<ArgsInfoType>::
+InitializeImageType()
+{
+ ADD_IMAGE_TYPE(Dim, uchar);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+void clitk::RelativePositionDataBaseBuilderGenericFilter<ArgsInfoType>::
+SetArgsInfo(const ArgsInfoType & a)
+{
+ mArgsInfo=a;
+ SetIOVerbose(mArgsInfo.verbose_flag);
+ if (mArgsInfo.imagetypes_flag) this->PrintAvailableImageTypes();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+// Update with the number of dimensions and the pixeltype
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class FilterType>
+void clitk::RelativePositionDataBaseBuilderGenericFilter<ArgsInfoType>::
+SetOptionsFromArgsInfoToFilter(FilterType * f)
+{
+ f->SetAFDBFilename(mArgsInfo.afdb_arg);
+ f->SetAFDBPath(mArgsInfo.afdb_path_arg);
+ f->SetNumberOfBins(mArgsInfo.bins_arg);
+ f->SetNumberOfAngles(mArgsInfo.nb_arg);
+ f->SetAreaLossTolerance(mArgsInfo.tol_arg);
+ f->SetSupportName(mArgsInfo.supportName_arg);
+ f->SetTargetName(mArgsInfo.targetName_arg);
+ for(int i=0; i<mArgsInfo.objectName_given; i++)
+ f->AddObjectName(mArgsInfo.objectName_arg[i]);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+// Update with the number of dimensions and the pixeltype
+//--------------------------------------------------------------------
+template<class ArgsInfoType>
+template<class ImageType>
+void clitk::RelativePositionDataBaseBuilderGenericFilter<ArgsInfoType>::
+UpdateWithInputImageType()
+{
+ // Create filter
+ typedef clitk::RelativePositionDataBaseBuilderFilter<ImageType> FilterType;
+ typename FilterType::Pointer filter = FilterType::New();
+
+ // Set global Options
+ SetOptionsFromArgsInfoToFilter<FilterType>(filter);
+
+ // Go !
+ filter->Update();
+
+ // Write/Save results
+ // typename ImageType::Pointer output = filter->GetOutput();
+ //this->template SetNextOutput<ImageType>(output);
+
+}
+//--------------------------------------------------------------------
+
+
- CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
===========================================================================**/
-#ifndef CLITKEXTRACTPATIENTGENERICFILTER_H
-#define CLITKEXTRACTPATIENTGENERICFILTER_H
+#ifndef CLITKRELATIVEPOSITIONGENERICFILTER_H
+#define CLITKRELATIVEPOSITIONGENERICFILTER_H
// clitk
#include "clitkIO.h"
#include "clitkRelativePositionGenericFilter.txx"
#endif
-#endif // #define CLITKEXTRACTPATIENTGENERICFILTER_H
+#endif // #define CLITKRELATIVEPOSITIONGENERICFILTER_H
f->SetVerboseStepFlag(mArgsInfo.verboseStep_flag);
f->SetWriteStepFlag(mArgsInfo.writeStep_flag);
+ // Must be set before AddOrientationTypeString
+ f->SetInverseOrientationFlag(mArgsInfo.inverse_flag);
+
for(uint i=0; i<mArgsInfo.orientation_given; i++) {
f->AddOrientationTypeString(mArgsInfo.orientation_arg[i]);
}
f->SetRemoveObjectFlag(!mArgsInfo.doNotRemoveObject_flag);
f->SetAutoCropFlag(!mArgsInfo.noAutoCrop_flag);
f->SetCombineWithOrFlag(mArgsInfo.combineWithOr_flag);
- f->SetInverseOrientationFlag(mArgsInfo.inverse_flag);
-
}
//--------------------------------------------------------------------
filter->SetInput(input);
filter->SetInputObject(object);
if (mArgsInfo.angle1_given && mArgsInfo.angle2_given)
- filter->AddAngles(mArgsInfo.angle1_arg, mArgsInfo.angle2_arg);
+ filter->AddAnglesInDeg(mArgsInfo.angle1_arg, mArgsInfo.angle2_arg);
SetOptionsFromArgsInfoToFilter<FilterType>(filter);
// Go !
FilterType::Pointer filter = FilterType::New();
filter->SetArgsInfo(args_info);
- filter->Update();
+ CLITK_TRY_CATCH_EXIT(filter->Update());
// this is the end my friend
return EXIT_SUCCESS;
clitk::SetBackgroundGenericFilter::Pointer genericFilter=clitk::SetBackgroundGenericFilter::New();
genericFilter->SetArgsInfo(args_info);
- genericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(genericFilter->Update());
return EXIT_SUCCESS;
}// end main
clitk::WarpImageGenericFilter::Pointer genericFilter=clitk::WarpImageGenericFilter::New();
genericFilter->SetArgsInfo(args_info);
- genericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(genericFilter->Update());
return EXIT_SUCCESS;
}// end main
zeroVFGenericFilter->SetVerbose(args_info.verbose_flag);
//update
- zeroVFGenericFilter->Update();
+ CLITK_TRY_CATCH_EXIT(zeroVFGenericFilter->Update());
return EXIT_SUCCESS;
}
#endif
git://git.creatis.insa-lyon.fr / clitk.git / commitdiff