# Extract labelmap contour

**URL:** <https://discourse.slicer.org/t/extract-labelmap-contour/410>\
**Category:** Development\
**Created:** [May 31, 2017, 10:09pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410 "2017-05-31T22:09:24Z")\
**Posts on this page:** 11\
**Page:** 1

<div class="post-metadata">

**Author:** ![DanC](https://avatars.discourse-cdn.com/v4/letter/d/77aa72/32.png) [@DanC](https://discourse.slicer.org/u/DanC)\
**Post date:** [May 31, 2017, 10:09pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/1 "2017-05-31T22:09:24Z")

</div>

Hi.  
I would like to extract the contour of the labelmap saved in a nrrd file using itk filters. It seems that the toggle button for the label map in the reslice viewers has this functionality. Is there a sample code of this functionality that I could search for? I’ve been trying to use the itk::BinaryContourImageFilter, but cannot get it to work.

Thanks.

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [June 2, 2017, 2:23pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/2 "2017-06-02T14:23:42Z")

</div>

Would you like an image, a closed surface, or a set of planar contours? In what file format?

---

<div class="post-metadata">

**Author:** ![DanC](https://avatars.discourse-cdn.com/v4/letter/d/77aa72/32.png) [@DanC](https://discourse.slicer.org/u/DanC)\
**Post date:** [June 2, 2017, 9:18pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/3 "2017-06-02T21:18:49Z")

</div>

I would like to extract the point set (cloud) of the label contours. Here is my current code.  
Thanks!

```
// declare the image class
typedef itk::Image<float, 3> OriginalImageType;

// -Read in label map
typedef itk::ImageFileReader<OriginalImageType> itkReaderType;
itkReaderType::Pointer itkreader = itkReaderType::New();
itkreader->SetFileName("label.nrrd");
itkreader->Update();

OriginalImageType::Pointer outImage = itkreader->GetOutput();

// Extract point data from label map
typedef itk::ImageRegionConstIterator< OriginalImageType > IteratorType;
typedef itk::PointSet< float, 3 > KdPointSetType;
typedef KdPointSetType::PointType KdPointType;

IteratorType it(outImage, outImage->GetBufferedRegion());
it.GoToBegin();

KdPointSetType::Pointer PointSet = KdPointSetType::New();
KdPointType point;

unsigned long pointId = 0;
	
while (!it.IsAtEnd())
{
	// Convert the pixel position into a Point
	outImage->TransformIndexToPhysicalPoint(it.GetIndex(), point);
	PointSet->SetPoint(pointId, point);
	// Transfer the pixel data to the value associated with the point.
	PointSet->SetPointData(pointId, it.Get());
	++it;
	++pointId;
}
std::cout << "Number Of Points = " << PointSet->GetNumberOfPoints() << std::endl;
```

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [June 2, 2017, 9:49pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/4 "2017-06-02T21:49:34Z")

</div>

If you only have a point cloud then you cannot do much with your data (because you don’t know what is inside/outside, you don’t have surface normals, you don’t have surfaces that you can render, etc). until you reconstruct a surface.

That’s why typically you run a marching cubes filter to create surface mesh from the image.

---

<div class="post-metadata">

**Author:** ![DanC](https://avatars.discourse-cdn.com/v4/letter/d/77aa72/32.png) [@DanC](https://discourse.slicer.org/u/DanC)\
**Post date:** [June 2, 2017, 10:01pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/5 "2017-06-02T22:01:27Z")

</div>

I am trying to write a deformablemesh surface filter that would deform against the point cloud of the surface. The slicers are so far apart that the marching cubes won’t give a satisfying result. Thanks for the help!

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [June 2, 2017, 10:20pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/6 "2017-06-02T22:20:26Z")

</div>

Reconstruction from point cloud is a much difficult problem than smooth interpolation between slices, so that’s not a good direction to explore.

Instead, create a high-resolution labelmap with isotropic spacing (preferably create a high-res labelmap in the first place by segmenting an isotropic supersampled input volume) using Crop volume module. For deformation, use Transforms module directly on the labelmap; or import the labelmap to a Segmentation node, create closed surface representation (point cloud with connectivity; you can apply smoothing on the labelmap in Segment editor; or during conversion to closed surface), make it the master representation and apply transform to that.

No need for any custom module development, these are all very common tasks.

---

<div class="post-metadata">

**Author:** ![Lorensen](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lorensen/32/31_2.png) [@Lorensen](https://discourse.slicer.org/u/Lorensen)\
**Post date:** [June 2, 2017, 11:42pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/7 "2017-06-02T23:42:24Z")

</div>

Vtk has recently added a number of point cloud filters. See  
[https://lorensen.github.io/VTKExamples/site/Cxx/#point-cloud-operations](https://lorensen.github.io/VTKExamples/site/Cxx/#point-cloud-operations)

In particular there is a new surface construction filter:  
[https://lorensen.github.io/VTKExamples/site/Cxx/Points/ExtractSurface/](https://lorensen.github.io/VTKExamples/site/Cxx/Points/ExtractSurface/)

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [June 3, 2017, 3:42am UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/8 "2017-06-03T03:42:28Z")

</div>

Thanks Bill for the links. These methods are useful when you don’t have connectivity information in your point cloud and you want to reconstruct only a single structure.

If you have connectivity information between the points (you have a surface mesh), then it would be a huge mistake to discard that information. Even when connectivity information within an image slice is available and you only miss connectivity info between slices (this is how segmentation is stored in DICOM-RT structure sets), there is a lot of ambiguity in reconstructing the surface. There is an [algorithm for this in SlicerRT](https://app.assembla.com/spaces/slicerrt/subversion/source/HEAD/trunk/SlicerRt/src/DicomRtImportExport/ConversionRules/vtkPlanarContourToClosedSurfaceConversionRule.cxx), that can reconstruct nice smooth surfaces most of the time while preserve most details, it’s quite tricky and in some cases the results are still not optimal.

I think any of these methods would give better results (no reconstruction errors, controllable smoothing, works well with multiple labels):

- Option A. Convert labelmap to mesh, smooth the mesh, then convert back to labelmap. This method is implemented in the Segment Editor: Smoothing effect, Joint smoothing method. This is somewhat similar to what @DanC may try to achieve, but with the important difference of not discarding connectivity information when we convert labelmap to mesh.
- Option B. Reduce the slice spacing to approximately match in-slice pixel size and insert empty slices between the original slices. Then run [morphological contour interpolation](https://github.com/KitwareMedical/ITKMorphologicalContourInterpolation). The morphological interpolator is available in Slicer in the Segment editor module, as `Interpolate between slices` effect.
- Option C. Supersample the labelmap and apply smoothing (median filter, morphological operators, Gaussian, etc). These are available in Segment editor’s smoothing effect.

---

<div class="post-metadata">

**Author:** ![DanC](https://avatars.discourse-cdn.com/v4/letter/d/77aa72/32.png) [@DanC](https://discourse.slicer.org/u/DanC)\
**Post date:** [June 3, 2017, 5:31pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/9 "2017-06-03T17:31:28Z")

</div>

Thanks for the suggestions Andras.

I am trying to reconstruct histological cross sections of a prostate including other structures such as urethra,etc.  
The spacing of the image set is 0.07 x 0.07 x 2 mm.  
The cross sections varies so much between images, that any smoothing or interpolation I’ve tried in slicer would not give a continuous surface.Hence I was trying to explore a totally different approach. The urethra would look like a spine instead of a continuous surface after reconstruction.

Option A. Are the suggested algorithms different than the ones used in the Model Maker Module? The smoothing/ joint smoothing in the model maker module wasn’t good enough in my case.  
Option B. I will try the morphological contour interpolation. I wan’t aware of this new function.  
Option C. Does supersample meen to resample with smaller spacing between slices with interpolation?

Thanks!

---

<div class="post-metadata">

**Author:** ![DanC](https://avatars.discourse-cdn.com/v4/letter/d/77aa72/32.png) [@DanC](https://discourse.slicer.org/u/DanC)\
**Post date:** [June 3, 2017, 5:33pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/10 "2017-06-03T17:33:56Z")

</div>

Hi Bill.  
Thank for the suggestion. I will definitely explore these if I go this direction.

Any suggestions to extract the contour points from the Slicer label maps? I currently am able to save the point set of all filled voxels using ITK.

Thanks!

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [June 3, 2017, 5:58pm UTC](https://discourse.slicer.org/t/extract-labelmap-contour/410/11 "2017-06-03T17:58:53Z")

</div>

> [@DanC](#):
>
> Are the suggested algorithms different than the ones used in the Model Maker Module?

The one in Segment Editor is slightly improved, but it still won’t be able to do any good if voxel spacing is highly anisotropic (you’ll see staircase artifacts or you’ll lose details in the reconstructed surface). **If spacing of your volume differs more than about 10% between axes then first you have to crop&resample your volume using [Crop Volume](https://www.slicer.org/wiki/Documentation/Nightly/Modules/CropVolume) module.**

So, before trying any of the recommended options, you must resample your input volume.

> [@DanC](#):
>
> Does supersample meen to resample with smaller spacing between slices with interpolation?

Yes.

When you have 20x difference between intra/inter-slice spacing, then most likely you cannot utilize all the information in the slice, so you have to use a spacing that is in between the smallest spacing vale - for 0.07x0.07x2mm spacing, something like 0.5mm should be good enough.

Do you register histology to US or MRI? Do you use the [SlicerProstate](https://www.slicer.org/wiki/Documentation/Nightly/Extensions/SlicerProstate) or [Segment registration](https://github.com/SlicerRt/SegmentRegistration) extensions? They are specifically developed for prostate registration between high-resolution and low-resolution images.