+++ /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
+++ /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
git://git.creatis.insa-lyon.fr / clitk.git / commitdiff