at.SetFromDataElement( contourdata );
const double* points = at.GetValues();
<<<<<<< Updated upstream
+<<<<<<< variant A
+>>>>>>> variant B
+ // unsigned int npts = at.GetNumberOfValues() / 3;
+======= end
=======
// unsigned int npts = at.GetNumberOfValues() / 3;
>>>>>>> Stashed changes
itkNewMacro(Self);
void Print(std::ostream & os = std::cout) const;
+
#if GDCM_MAJOR_VERSION == 2
- bool Read(gdcm::Item const & item);
+ bool Read(gdcm::Item * item);
+ void UpdateDicomItem();
#else
bool Read(gdcm::SQItem * item);
#endif
+
vtkPolyData * GetMesh();
+ void SetMesh(vtkPolyData * mesh);
vtkPoints * GetPoints() {return mData;}
double GetZ() const {return mZ;}
+
protected:
- void ComputeMesh();
+ void ComputeMeshFromDataPoints();
+ void ComputeDataPointsFromMesh();
unsigned int mNbOfPoints;
std::string mType;
vtkSmartPointer<vtkPoints> mData;
bool mMeshIsUpToDate;
///Z location of the contour
double mZ;
+
+ gdcm::Item * mItem;
private:
DicomRT_Contour();
#include "clitkDicomRT_ROI.h"
#include <vtkSmartPointer.h>
#include <vtkAppendPolyData.h>
+#include <vtkImageClip.h>
+#include <vtkMarchingSquares.h>
#if GDCM_MAJOR_VERSION == 2
#include "gdcmAttribute.h"
{
mName = "NoName";
mNumber = -1;
+ mImage = NULL;
mColor.resize(3);
mColor[0] = mColor[1] = mColor[2] = 0;
mMeshIsUpToDate = false;
return false;
}
const gdcm::DataElement& csq = nestedds.GetDataElement( tcsq );
- gdcm::SmartPointer<gdcm::SequenceOfItems> sqi2 = csq.GetValueAsSQ();
+ mContoursSequenceOfItems = csq.GetValueAsSQ();
+ gdcm::SmartPointer<gdcm::SequenceOfItems> & sqi2 = mContoursSequenceOfItems;
if( !sqi2 || !sqi2->GetNumberOfItems() )
{
}
for(unsigned int i = 0; i < nitems; ++i)
{
- const gdcm::Item & j = sqi2->GetItem(i+1); // Item start at #1
+ gdcm::Item & j = sqi2->GetItem(i+1); // Item start at #1
DicomRT_Contour::Pointer c = DicomRT_Contour::New();
- bool b = c->Read(j);
+ bool b = c->Read(&j);
if (b) {
mListOfContours.push_back(c);
}
vtkPolyData * clitk::DicomRT_ROI::GetMesh()
{
if (!mMeshIsUpToDate) {
- ComputeMesh();
+ ComputeMeshFromContour();
}
return mMesh;
}
//--------------------------------------------------------------------
-void clitk::DicomRT_ROI::ComputeMesh()
+void clitk::DicomRT_ROI::ComputeMeshFromContour()
{
vtkSmartPointer<vtkAppendPolyData> append = vtkSmartPointer<vtkAppendPolyData>::New();
for(unsigned int i=0; i<mListOfContours.size(); i++) {
gdcm::DataSet & ds = mItemInfo->GetNestedDataSet();
ds.Replace(de);
+ // From MESH to CONTOURS
+ ComputeContoursFromImage();
+
// Update contours
+ DD(mListOfContours.size());
for(uint i=0; i<mListOfContours.size(); i++) {
DD(i);
- DicomRT_Contour::Pointer roi = mListOfContours[i];
- roi->UpdateDicomItem(mItemContour);
+ DicomRT_Contour::Pointer contour = mListOfContours[i];
+ contour->UpdateDicomItem();//mItemContour);
}
+
+ // Nb of contours
+ unsigned int nitems = mContoursSequenceOfItems->GetNumberOfItems();
+ DD(nitems);
+
+ // Write [Contour Sequence] = 0x3006,0x0040)
+ gdcm::DataSet & dsc = mItemContour->GetNestedDataSet();
+ gdcm::Tag tcsq(0x3006,0x0040);
+ const gdcm::DataElement& csq = dsc.GetDataElement( tcsq );
+ gdcm::DataElement dec(csq);
+ dec.SetValue(*mContoursSequenceOfItems);
+ dsc.Replace(dec);
+
+ gdcm::DataSet & a = mContoursSequenceOfItems->GetItem(1).GetNestedDataSet();
+ gdcm::Attribute<0x3006,0x0050> at;
+ gdcm::Tag tcontourdata(0x3006,0x0050);
+ gdcm::DataElement contourdata = a.GetDataElement( tcontourdata );
+ at.SetFromDataElement( contourdata );
+ const double* points = at.GetValues();
+ DD(points[0]);
+
}
//--------------------------------------------------------------------
#endif
return mImage;
}
//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+void clitk::DicomRT_ROI::ComputeContoursFromImage()
+{
+ DD("ComputeMeshFromImage");
+
+ // Check that an image is loaded
+ if (!mImage) return;
+
+ // Only consider 3D here
+ if (mImage->GetNumberOfDimensions() != 3) {
+ FATAL("DicomRT_ROI::ComputeMeshFromImage only work with 3D images");
+ }
+
+ // Get the VTK image
+ vtkImageData * image = mImage->GetVTKImages()[0];
+
+ // Get initial extend for the clipping
+ vtkSmartPointer<vtkImageClip> clipper = vtkSmartPointer<vtkImageClip>::New();
+ clipper->SetInput(image);
+ int* extent = image->GetExtent();
+ DDV(extent, 6);
+ // std::vector<int> extend;
+
+ // Prepare the marching squares
+ vtkSmartPointer<vtkMarchingSquares> squares = vtkSmartPointer<vtkMarchingSquares>::New();
+ squares->SetInput(clipper->GetOutput());
+
+ // Loop on slice
+ uint n = image->GetDimensions()[2];
+ DD(n);
+ for(uint i=0; i<n; i++) {
+ // Clip to the current slice
+ extent[4] = extent[5] = image->GetOrigin()[2]+i*image->GetSpacing()[2];
+ DDV(extent, 6);
+ clipper->SetOutputWholeExtent(extent[0],extent[1],extent[2],
+ extent[3],extent[4],extent[5]);
+
+
+ squares->Update();
+ DD(squares->GetNumberOfContours());
+ mListOfContours[i]->SetMesh(squares->GetOutput());
+ }
+}
+//--------------------------------------------------------------------
DicomRT_Contour* GetContour(int n);
// Compute a vtk mesh from the dicom contours
- void ComputeMesh();
+ void ComputeMeshFromContour();
+ void ComputeContoursFromImage();
// Indicate if the mesh is uptodate according to the dicom
void SetDicomUptodateFlag(bool b) { m_DicomUptodateFlag = b; }
#if GDCM_MAJOR_VERSION == 2
gdcm::Item * mItemInfo;
gdcm::Item * mItemContour;
+ gdcm::SmartPointer<gdcm::SequenceOfItems> mContoursSequenceOfItems;
#endif
private:
const gdcm::DataElement &roicsq = ds.GetDataElement( troicsq );
gdcm::DataElement de2(roicsq);
de2.SetValue(*mROIContoursSequenceOfItems);
- ds.Replace(de);
-
+ ds.Replace(de2);
+
+ //DEBUG
+ gdcm::DataSet & a = mROIContoursSequenceOfItems->GetItem(1).GetNestedDataSet();
+ gdcm::Tag tcsq(0x3006,0x0040);
+ const gdcm::DataElement& csq = a.GetDataElement( tcsq );
+ gdcm::SmartPointer<gdcm::SequenceOfItems> sqi2 = csq.GetValueAsSQ();
+ gdcm::Item & j = sqi2->GetItem(1);
+ gdcm::DataSet & b = j.GetNestedDataSet();
+ gdcm::Attribute<0x3006,0x0050> at;
+ gdcm::Tag tcontourdata(0x3006,0x0050);
+ gdcm::DataElement contourdata = b.GetDataElement( tcontourdata );
+ at.SetFromDataElement( contourdata );
+ const double* points = at.GetValues();
+ DD(points[0]);
+
+
// Write dicom
gdcm::Writer writer;
//writer.CheckFileMetaInformationOff();
// clitk
#include "clitkImage2DicomRTStructFilter.h"
#include "clitkImageCommon.h"
+#include "vvFromITK.h"
// vtk
#include <vtkPolyDataToImageStencil.h>
#include <vtkLinearExtrusionFilter.h>
#include <vtkMetaImageWriter.h>
#include <vtkImageData.h>
+
+// itk
#include <itkImage.h>
#include <itkVTKImageToImageFilter.h>
template<class PixelType>
void clitk::Image2DicomRTStructFilter<PixelType>::Update()
{
- DD("Image2DicomRTStructFilter::Update");
-
-
+ DD("Image2DicomRTStructFilter::GenerateData");
+
+ // Read DicomRTStruct
+ std::string filename = "RS.zzQAnotmt_french01_.dcm";
+ clitk::DicomRT_StructureSet::Pointer structset = clitk::DicomRT_StructureSet::New();
+ structset->Read(filename);
+
+ DD(structset->GetName());
+ clitk::DicomRT_ROI * roi = structset->GetROIFromROINumber(1); // Aorta
+ DD(roi->GetName());
+ DD(roi->GetROINumber());
+
+ // Add an image to the roi
+ vvImage::Pointer im = vvImageFromITK<3, PixelType>(m_Input);
+ roi->SetImage(im);
+
+ // Get one contour
+ DD("Compute Mesh");
+ roi->ComputeMeshFromImage();
+ vtkSmartPointer<vtkPolyData> mesh = roi->GetMesh();
+ DD("done");
+
+ // Change the mesh (shift by 10);
+ // const vtkSmartPointer<vtkPoints> & points = mesh->GetPoints();
+ // for(uint i=0; i<mesh->GetNumberOfVerts (); i++) {
+ // DD(i);
+ // double * p = points->GetPoint(i);
+ // p[0] += 30;
+ // points->SetPoint(i, p);
+ // }
+ roi->SetName("TOTO");
+ roi->SetDicomUptodateFlag(false); // indicate that dicom info must be updated from the mesh.
+
+ // Convert to dicom ?
+ DD("TODO");
+
+ // Write
+ structset->Write("toto.dcm");
}
//--------------------------------------------------------------------
ExtractStation_3A()
{
if (!CheckForStation("3A")) return;
-
+
+ StartNewStep("Station 3A");
+ StartSubStep();
+
// Get the current support
StartNewStep("[Station 3A] Get the current 3A suppport");
m_Working_Support = m_ListOfSupports["S3A"];
ExtractStation_3A_AntPost_S5();
ExtractStation_3A_AntPost_S6();
ExtractStation_3A_AntPost_Superiorly();
-
ExtractStation_3A_Remove_Structures();
Remove_Structures("3A", "Aorta");
Remove_Structures("3A", "CommonCarotidArteryLeft");
Remove_Structures("3A", "CommonCarotidArteryRight");
Remove_Structures("3A", "BrachioCephalicArtery");
- // Remove_Structures("3A", "BrachioCephalicVein"); ?
+ ExtractStation_3A_PostToBones();
+
+
// ExtractStation_3A_Ant_Limits(); --> No, already in support; to remove
// ExtractStation_3A_Post_Limits(); --> No, more complex, see Vessels etc
writeImage<MaskImageType>(m_ListOfStations["3A"], "seg/Station3A.mhd");
GetAFDB()->SetImageFilename("Station3A", "seg/Station3A.mhd");
WriteAFDB();
+ StopSubStep();
}
//--------------------------------------------------------------------
MaskImagePointer SVC = GetAFDB()->template GetImage <MaskImageType>("SVC");
// Trial in 3D -> difficulties superiorly. Stay slice by slice.
- /*
- typedef clitk::AddRelativePositionConstraintToLabelImageFilter<MaskImageType> RelPosFilterType;
- typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New();
- relPosFilter->VerboseStepFlagOff();
- relPosFilter->WriteStepFlagOff();
- relPosFilter->SetBackgroundValue(GetBackgroundValue());
- relPosFilter->SetInput(m_Working_Support);
- relPosFilter->SetInputObject(SVC);
- relPosFilter->AddOrientationTypeString("PostTo");
- relPosFilter->SetInverseOrientationFlag(true);
- relPosFilter->SetIntermediateSpacing(4);
- relPosFilter->SetIntermediateSpacingFlag(false);
- relPosFilter->SetFuzzyThreshold(0.5);
- // relPosFilter->AutoCropFlagOff(); // important ! because we join the slices after this loop
- relPosFilter->Update();
- m_Working_Support = relPosFilter->GetOutput();
- */
-
// Slice by slice not post to SVC. Use initial spacing
m_Working_Support =
clitk::SliceBySliceRelativePosition<MaskImageType>(m_Working_Support, SVC, 2,
// resize like support, extract slices
// while single CCL -> remove
// when two remove only the most post
+ MaskImagePointer BrachioCephalicVein =
+ GetAFDB()->template GetImage <MaskImageType>("BrachioCephalicVein");
BrachioCephalicVein = clitk::ResizeImageLike<MaskImageType>(BrachioCephalicVein,
m_Working_Support,
GetBackgroundValue());
for(uint i=0; i<slices.size(); i++) {
// Labelize slices_BCV
slices_BCV[i] = Labelize<MaskSliceType>(slices_BCV[i], 0, true, 1);
+
// Compute centroids
std::vector<typename MaskSliceType::PointType> centroids;
ComputeCentroids<MaskSliceType>(slices_BCV[i], GetBackgroundValue(), centroids);
- if (centroids.size() > 1) {
+
+ // If several centroid, select the one most anterior
+ if (centroids.size() > 2) {
// Only keep the one most post
typename MaskSliceType::PixelType label;
if (centroids[1][1] > centroids[2][1]) {
- label = 1;
- }
- else {
label = 2;
}
-
- HERE
-
- slices_BCV[i] = clitk::SetBackground<MaskSliceType>(slices_BCV[i], slices_BCV[i],
- label,
- GetBackgroundValue(), true);
+ else {
+ label = 1;
+ }
+ // "remove" the CCL
+ slices_BCV[i] = clitk::SetBackground<MaskSliceType, MaskSliceType>(slices_BCV[i],
+ slices_BCV[i],
+ label,
+ GetBackgroundValue(),
+ true);
}
+
// Remove from the support
- slices[i] = clitk::AndNot(slices[i], slices_BCV[i], GetBackgroundValue());
+ clitk::AndNot<MaskSliceType>(slices[i], slices_BCV[i], GetBackgroundValue());
}
// Joint
{
if (!CheckForStation("3P")) return;
+ StartNewStep("Station 3P");
+ StartSubStep();
+
// Get the current support
StartNewStep("[Station 3P] Get the current 3P suppport");
m_Working_Support = m_ListOfSupports["S3P"];
writeImage<MaskImageType>(m_ListOfStations["3P"], "seg/Station3P.mhd");
GetAFDB()->SetImageFilename("Station3P", "seg/Station3P.mhd");
WriteAFDB();
+ StopSubStep();
}
//--------------------------------------------------------------------
clitk::ExtractLymphStationsFilter<TImageType>::
ExtractStation_8()
{
- if (CheckForStation("8")) {
- ExtractStation_8_SI_Limits(); // OK, validated
- ExtractStation_8_Ant_Limits(); // OK, validated
- ExtractStation_8_Left_Sup_Limits(); // OK, validated
- ExtractStation_8_Left_Inf_Limits(); // OK, validated
- ExtractStation_8_Single_CCL_Limits(); // OK, validated
- ExtractStation_8_Remove_Structures(); // OK, validated
-
- // Store image filenames into AFDB
- writeImage<MaskImageType>(m_ListOfStations["8"], "seg/Station8.mhd");
- GetAFDB()->SetImageFilename("Station8", "seg/Station8.mhd");
- WriteAFDB();
- }
+ if (!CheckForStation("8")) return;
+
+ StartNewStep("Station 8");
+ StartSubStep();
+ ExtractStation_8_SI_Limits(); // OK, validated
+ ExtractStation_8_Ant_Limits(); // OK, validated
+ ExtractStation_8_Left_Sup_Limits(); // OK, validated
+ ExtractStation_8_Left_Inf_Limits(); // OK, validated
+ ExtractStation_8_Single_CCL_Limits(); // OK, validated
+ ExtractStation_8_Remove_Structures(); // OK, validated
+
+ // Store image filenames into AFDB
+ writeImage<MaskImageType>(m_ListOfStations["8"], "seg/Station8.mhd");
+ GetAFDB()->SetImageFilename("Station8", "seg/Station8.mhd");
+ WriteAFDB();
+ StopSubStep();
}
//--------------------------------------------------------------------
void ExtractStation_3A_AntPost_S5();
void ExtractStation_3A_AntPost_S6();
void ExtractStation_3A_AntPost_Superiorly();
+ void ExtractStation_3A_Remove_Structures();
+
void ExtractStation_3A_SI_Limits();
void ExtractStation_3A_Ant_Limits();
void ExtractStation_3A_Post_Limits();
}
// Extract Station8
- StartNewStep("Station 8");
- StartSubStep();
ExtractStation_8();
- StopSubStep();
// Extract Station3P
- StartNewStep("Station 3P");
- StartSubStep();
ExtractStation_3P();
- StopSubStep();
// Extract Station3A
- StartNewStep("Station 3A");
- StartSubStep();
ExtractStation_3A();
- StopSubStep();
-
// HERE
// Write result
clitk::DicomRT_StructureSet::Pointer s = filter.GetDicomRTStruct();
- s->Write(args_info.output_arg);
+ // s->Write(args_info.output_arg);
// This is the end my friend
return 0;
}
+//--------------------------------------------------------------------