VTK/Closed Surface Clipping: Difference between revisions

From KitwarePublic
< VTK
Jump to navigationJump to search
 
(69 intermediate revisions by the same user not shown)
Line 1: Line 1:
__NOTOC__
__NOTOC__
[[Image:ClipClosedSurface.png|right|300px|ClipClosedSurface regression test.]]
[[Image:ClipClosedSurface.png|right|300px|thumb|A clipped surface with generation of color scalars.]]
The main goal of this project is to make it possible to clip a surface, and then create a "cap" polygon (or polygons) that fill in the area that is enclosed by all the new edges that were created by the clipping.  In other words, the goal is to close the surface after the clipping has been performed.  For this purpose, a new filter called vtkClipClosedSurface has been created.
The main goal of this project is to make it possible to clip a surface, and then create a "cap" polygon (or polygons) that fill in the area that is enclosed by all the new edges that were created by the clipping.  In other words, the goal is to close the surface after the clipping has been performed.  For this purpose, a new filter called vtkClipClosedSurface has been created.


This class was contributed to VTK by David Gobbi on March 30, 2010.
This class was committed to VTK cvs by David Gobbi on March 30, 2010. It is in VTK 5.8 and later releases.
 
== Features ==
== Features ==
 
[[Image:ClipSubtractedSurface.png|right|300px|thumb|A skin isosurface with a reversed bone isosurface.]]
The vtkClipClosedSurface filter has the following features:
The vtkClipClosedSurface filter has the following features:
* A vtkPlaneCollection is used to set the clip region (arbitrary implicit functions are not allowed).
* A vtkPlaneCollection is used to set the clip region (arbitrary implicit functions are not allowed).
* If given a closed surface as input, it will produce a closed surface as output.
* When given a closed surface as input, it will produce a closed surface as output.
* If given an open surface as input, it will create a closed surface if supplied with a clipping plane that clips away all open edges.
* When given an open surface as input, it will create a closed surface if all of the boundary edges lie outside of the clip region.
* The input surface can be concave, and can even be composed of multiple closed surfaces, as long as there is no intersection.
* The input surface can be concave, but cannot self-intersect.
* The input can have multiple surfaces, where any surfaces that are inside-out will be treated as holes within larger, enclosing surfaces.  The inner and outer surfaces must not intersect each other.
* If an input surface is inside-out, but does not have a larger surface enclosing it, then the cut faces will be left open.
 
Note: an inside-out surface can be created by applying vtkReverseSense to a regular surface.


== Caveats ==
== Caveats ==
The input surface must satisfy these conditions:
* It should not self-intersect (polygons can only touch at their edges).
* It should not have open edges, but the filter will succeed if all open edges are clipped away.
* It should not have any non-manifold edges (edges that are shared by more than two faces).
* (Note: vtkFeatureEdges can be used to test the edges of a data set).


The newly-created clip surfaces have the following limitations:
The newly-created clip surfaces have the following limitations:
* If the contours at the clip plane are not closed, or if they don't have simple topology, an incorrect surface will be generated.
* The triangles that are generated are often very narrow.  The code does not attempt to generate high-quality triangles.
* The triangles that are generated are often very narrow.  The code does not attempt to generate high-quality triangles.
* This filter ignores vtkVertex cells, and the output will contain no verts.
* This filter ignores vtkVertex cells, and the output will contain no verts.
Line 23: Line 31:
== Options ==
== Options ==


[[Image:ClippedVolumeOutline.png|right|300px|ClippedVolumeOutline]]
[[Image:ClippedVolumeOutline.png|right|300px|thumb|Clipping a cube and displaying the outline only.]]
* '''GenerateFaces''' - On by default.  If turned off, the output will have no polygons.  This is useful if only the outlines are desired.
* '''Tolerance''' - 1e-6 by default.  This is the tolerance for merging points while clipping.
* '''TriangulationErrorDisplay''' - Off by default.  This will cause errors to be generated if triangulation fails.  Triangulation failures cause tiny holes and overlaps in the surface.
* '''GenerateFaces''' - On by default.  If turned off, the output will have no polygons, which useful if only the outlines are desired.
* '''GenerateOutline''' - Off by default.  If turned on, contours will be generated where the planes intersect the surface.
* '''GenerateOutline''' - Off by default.  If turned on, contours will be generated where the planes intersect the surface.
* '''GenerateColorScalars''' - Off by default.  Generate color scalars for the output, see below for more.
* '''ScalarMode''' - "None" by default, but can be set to "Colors" or "Labels".  Generate scalars for the output, see below for more.
* '''PassPointData''' - Off by default.  Turn on in order to pass point scalars, normals, and other point data to the output.


The '''GenerateColorScalars''' option creates color cell scalars so that the cut surfaces are displayed with a different color than the rest of the surface.  Cell scalars are used instead of point scalars because the use of point scalars would require point duplication where the cut surfaces meet the rest of the surface.  The lines cannot be colored differently from polys.  If you need the lines to be a different color, then you could use two instances of vtkClipClosedSurface: one to generates faces, and another to generates outlines.
Using '''SetScalarModeToColors''' or '''SetScalarModeToLabels''' will cause cell scalars to be generated so that the cut surfaces are displayed with a different color than the rest of the surface.  Cell scalars are used instead of point scalars because the use of point scalars would require point duplication where the cut surfaces meet the rest of the surface.  If GenerateOutline is on, then the lines will get the same scalar values as the faces.


If GenerateColorScalars is on, then three colors can be set:
If the ScalarMode is '''Labels''', then the original cells will have a scalar value of 0, cells on the clip surfaces will have a scalar value of 1, and cells on the plane with the ActivePlaneId will have a scalar value of 2.
 
If the ScalarMode is '''Colors''', then three colors can be set:
* '''ClipColor''' for the clip surfaces
* '''ClipColor''' for the clip surfaces
* '''BaseColor''' for the rest of the surface
* '''BaseColor''' for the rest of the surface
Line 37: Line 50:
The purpose of '''ActivePlaneColor''' is to color one cut plane differently, e.g. the plane that the user is interacting with.  The '''ActivePlaneId''' can be set to specify which of the cut planes is the active plane.
The purpose of '''ActivePlaneColor''' is to color one cut plane differently, e.g. the plane that the user is interacting with.  The '''ActivePlaneId''' can be set to specify which of the cut planes is the active plane.


The '''BaseColor''' will be ignored if the input already has 3-component color scalars.  If this is the case, the original scalars will be used instead of BaseColor.  This is useful if you want to chain two vtkClipClosedSurface filters together, or if you want to use vtkClipClosedSurface on the output of vtkOutlineSource or vtkVolumeOutlineSource.
The '''BaseColor''' will be ignored if the input already has 3-component color cell scalars.  If this is the case, the original scalars will be used instead of BaseColor.  This is useful if you want to chain two vtkClipClosedSurface filters together, or if you want to use vtkClipClosedSurface on the output of vtkOutlineSource or vtkVolumeOutlineSource.


Note that many VTK filters that create new cells (e.g. vtkStripper) will erase the cell scalars.  So if you use GenerateCellScalars, you should avoid using vtkStripper.
Note that many VTK filters that create new cells (e.g. vtkStripper) will erase the cell scalars.  So if you set ScalarMode, you should avoid using vtkStripper to post-process the data.


== Algorithm ==
== Algorithm ==


[[Image:ClipClosedSurfaceWireframe.png|right|300px|ClipClosedSurface wireframe.]]
[[Image:ClipClosedSurfaceWireframe.png|right|300px|thumb|Wireframe rendering, to show triangle quality.]]
This filter uses many of the functions that are built into vtkCell and its subclasses, but it also has several internal functions for polygon manipulation.  The tolerance for merging points etc. is set at d*1e-5 where "d" is the maximum dimension of the input bounds.  An overall description of what is done for each of the clipping planes is as follows:
This filter uses many of the functions that are built into vtkCell and its subclasses, but it also has several internal functions for polygon manipulation.  An overall description of what is done for each of the clipping planes is as follows:
# Clip scalars are generated for all points.
# Clip scalars are generated for all points.
# The data is clipped and contoured, using the clip scalars.
# The data is clipped and contoured, using the clip scalars.
Line 52: Line 65:
Generation of the cut faces from the cut contours is a multi-step process:
Generation of the cut faces from the cut contours is a multi-step process:
# The contour lines are joined end-to-end to form polygons.
# The contour lines are joined end-to-end to form polygons.
# If the polygons self-intersect, they are split at the intersection points.
# Any polygon points that form a 180 degree vertex are removed.
# Any polygon points that aren't at a corner of the polygon are removed.
# The sense of the polygons is checked against the cut plane normal to find out which polygons are holes.
# The sense of the polygons is set according to the cut plane normal.
# Two cuts are made between each hole and its containing polygon in order to generate two simple polygons.
# If any polygons are within other polygons, cuts are made between the outer and inner polygons.
# The polygons are checked for points where they pinch in on themselves, and are split at such points.
# The ear-cut triangulation method of vtkPolygon is used to triangulate the polygons.
# An ear-cut triangulation algorithm is used to triangulate the polygons.
# Some of the new triangles are subdivided in order to add back the points that were removed in step 3.
# Some of the new triangles are subdivided in order to add back the points that were removed in step 2.
 
== Future Work ==
 
=== Refactoring into related classes ===
 
==== vtkContourTriangulator (done, part of VTK 6.0) ====
 
This class would take a collection of contours (all within the same plane) and generate an appropriate set of polygons.  Around half of the code in vtkClipClosedSurface could be put into this class with very little modification.  This class would require the contours to be properly connected, so people would be advised to pass their line segments through vtkTriangleFilter before applying it.
 
This class would be used to do similar things as a constrained vtkDelaunay2D, but would not optimize the triangulation and would not be as sensitive to numerical instability.
 
==== vtkCutClosedSurface ====
 
This class could be responsible for generating just the cap polygons.  On the other hand, this could be provided as an option of the current class.
 
=== Guaranteeing a watertight output ===
 
==== Problem ====
 
The vtkFeatureEdges filter can check a polydata to see if it is watertight: a watertight surface will have no boundary edges and no non-manifold edges.  The vtkClipClosedSurface filter should guarantee in the strictest sense that it will produce a watertight surface if given a watertight surface as an input.  Unfortunately, VTK's low-level clipping and triangulation code makes this difficult to achieve.
 
# VTK clipping, as implemented in the vtkCell subclasses, requires the use of a point locator (usually vtkMergePoints) that merges coincident points.  Unfortunately, point merging can change connectivity and therefore change topology.  Specifically, a manifold data set can, after point merging, become non-manifold.
# VTK triangulation, as implemented in vtkPolygon, will leave holes rather than create degenerate or nearly-degenerate triangles.  So if the clipping produces polygons that (a) have any 180 degree vertices or (b) have vertices that are very close to each other, then the output will have free edges.
 
==== Solution ====
 
Unfortunately, the only solution is to do both of these:
# Write new clipping code that is specific to this filter, instead of using vtkCell clipping (DONE)
# Write new triangulation code, instead of using vtkPolygon::Triangulate()
 
Another reason not to use vtkPolygon::Triangulate() is performance: vtkPolygon's ear-cut algorithm has a worst-case performance of O(n^2), it was never meant for use on large polygons.  A simple locator-based ear-cut algorithm would be O(n) at best and O(n*log n) at worst.
 
A partial solution is to use a tolerance for creating new points along clipped edges. When an edge is interpolated during the clipping process, if the new point is within tolerance of either of the edge end points, then that endpoint will be used instead of the new point (DONE).  Also, if the clipped polygon is found to be degenerate, then the entire polygon can be removed at this stage so that it does not cause problems later.
 
=== Reducing the number of triangles produced ===
 
==== Problem ====
 
This class uses the clipping/contouring methods in the vtkCell subclasses, which triangulate polys while clipping or contouring them.  This produces many more points and cells than are actually needed to describe the geometry.  Also, when multiple cut planes are used, these cells can be clipped multiple times, causing a proliferation of small, sharp triangles.
 
==== Solution ====
 
There are two parts to the solution:
# Write new clipping/contouring code that doesn't triangulate the polys while clipping and contouring them. (DONE)
# Improve the triangulation of the cap polygons, or even go as far as not triangulating them at all.
 
The ideal solution is to write clipping/contouring code that can deal with concave polygons.  Then, in RequestData where the filter applies each clipping plane in sequence, the code can deal directly with the concave polygons produced by the previous clip-plane, and triangulation of the polygons can be deferred to the very end of RequestData.  The result will be that the final output will have the minimum possible number of triangles.
 
There is another good reason for writing specialized clipping/contouring code: it allows the orientation of the contours to be correctly set from the orientation of the clipped polygons.  This will eliminate the need to do an orientation check while creating the cap polygons, and will further guarantee that correct cap polygons are created.  For instance, right now vtkClipClosedSurface will create incorrectly-oriented cap polygons if its input is an inside-out surface.  But if the orientation of the cap polygons is set from the orientation of the original polygons (instead of being set according to the clip-plane normal like they are now), then cap polygons will be correctly generated even for inside-out surfaces. (DONE)
 
=== Handling non-closed inputs ===
 
==== Problem ====
 
The algorithm for creating new polygons assumes that all contours are closed.  Currently, when a non-closed contour is found, it's two ends are connected in order to close it.  This produces some odd results when the filter is given non-closed inputs.
 
==== Solution ====
 
There are three possible solutions:
# Generate an error when non-closed contours are found, but otherwise continue with the current behaviour.
# Throw out these non-closed contours.  Probably not a good idea.
# Add heuristics for optimally joining non-closed contours. (DONE)
 
The following heuristics have been added for joining the contours:
# When joining contours, the new edge must be on the hull of the total point set of all contours.
# The new edge must be of the correct sense, i.e. counterclockwise around the hull.
# If there are multiple candidate new edges for a loose end, the shortest is chosen.
# If there are no new edge candidates for a loose end, then that contour is removed.
 
Only joins along the hull are permitted, since it is assumed that the loose ends are the result of contours that intersect the bounds of the data set.

Latest revision as of 17:39, 1 February 2016

A clipped surface with generation of color scalars.

The main goal of this project is to make it possible to clip a surface, and then create a "cap" polygon (or polygons) that fill in the area that is enclosed by all the new edges that were created by the clipping. In other words, the goal is to close the surface after the clipping has been performed. For this purpose, a new filter called vtkClipClosedSurface has been created.

This class was committed to VTK cvs by David Gobbi on March 30, 2010. It is in VTK 5.8 and later releases.

Features

A skin isosurface with a reversed bone isosurface.

The vtkClipClosedSurface filter has the following features:

  • A vtkPlaneCollection is used to set the clip region (arbitrary implicit functions are not allowed).
  • When given a closed surface as input, it will produce a closed surface as output.
  • When given an open surface as input, it will create a closed surface if all of the boundary edges lie outside of the clip region.
  • The input surface can be concave, but cannot self-intersect.
  • The input can have multiple surfaces, where any surfaces that are inside-out will be treated as holes within larger, enclosing surfaces. The inner and outer surfaces must not intersect each other.
  • If an input surface is inside-out, but does not have a larger surface enclosing it, then the cut faces will be left open.

Note: an inside-out surface can be created by applying vtkReverseSense to a regular surface.

Caveats

The input surface must satisfy these conditions:

  • It should not self-intersect (polygons can only touch at their edges).
  • It should not have open edges, but the filter will succeed if all open edges are clipped away.
  • It should not have any non-manifold edges (edges that are shared by more than two faces).
  • (Note: vtkFeatureEdges can be used to test the edges of a data set).

The newly-created clip surfaces have the following limitations:

  • The triangles that are generated are often very narrow. The code does not attempt to generate high-quality triangles.
  • This filter ignores vtkVertex cells, and the output will contain no verts.
  • Similarly, input lines will be discarded unless the GenerateOutline option is on.

Options

Clipping a cube and displaying the outline only.
  • Tolerance - 1e-6 by default. This is the tolerance for merging points while clipping.
  • TriangulationErrorDisplay - Off by default. This will cause errors to be generated if triangulation fails. Triangulation failures cause tiny holes and overlaps in the surface.
  • GenerateFaces - On by default. If turned off, the output will have no polygons, which useful if only the outlines are desired.
  • GenerateOutline - Off by default. If turned on, contours will be generated where the planes intersect the surface.
  • ScalarMode - "None" by default, but can be set to "Colors" or "Labels". Generate scalars for the output, see below for more.
  • PassPointData - Off by default. Turn on in order to pass point scalars, normals, and other point data to the output.

Using SetScalarModeToColors or SetScalarModeToLabels will cause cell scalars to be generated so that the cut surfaces are displayed with a different color than the rest of the surface. Cell scalars are used instead of point scalars because the use of point scalars would require point duplication where the cut surfaces meet the rest of the surface. If GenerateOutline is on, then the lines will get the same scalar values as the faces.

If the ScalarMode is Labels, then the original cells will have a scalar value of 0, cells on the clip surfaces will have a scalar value of 1, and cells on the plane with the ActivePlaneId will have a scalar value of 2.

If the ScalarMode is Colors, then three colors can be set:

  • ClipColor for the clip surfaces
  • BaseColor for the rest of the surface
  • ActivePlaneColor to color one clip surface in a different color

The purpose of ActivePlaneColor is to color one cut plane differently, e.g. the plane that the user is interacting with. The ActivePlaneId can be set to specify which of the cut planes is the active plane.

The BaseColor will be ignored if the input already has 3-component color cell scalars. If this is the case, the original scalars will be used instead of BaseColor. This is useful if you want to chain two vtkClipClosedSurface filters together, or if you want to use vtkClipClosedSurface on the output of vtkOutlineSource or vtkVolumeOutlineSource.

Note that many VTK filters that create new cells (e.g. vtkStripper) will erase the cell scalars. So if you set ScalarMode, you should avoid using vtkStripper to post-process the data.

Algorithm

Wireframe rendering, to show triangle quality.

This filter uses many of the functions that are built into vtkCell and its subclasses, but it also has several internal functions for polygon manipulation. An overall description of what is done for each of the clipping planes is as follows:

  1. Clip scalars are generated for all points.
  2. The data is clipped and contoured, using the clip scalars.
  3. The contours are used to generate the cut faces.
  4. The process is repeated for the next cut plane

Generation of the cut faces from the cut contours is a multi-step process:

  1. The contour lines are joined end-to-end to form polygons.
  2. Any polygon points that form a 180 degree vertex are removed.
  3. The sense of the polygons is checked against the cut plane normal to find out which polygons are holes.
  4. Two cuts are made between each hole and its containing polygon in order to generate two simple polygons.
  5. The polygons are checked for points where they pinch in on themselves, and are split at such points.
  6. An ear-cut triangulation algorithm is used to triangulate the polygons.
  7. Some of the new triangles are subdivided in order to add back the points that were removed in step 2.

Future Work

Refactoring into related classes

vtkContourTriangulator (done, part of VTK 6.0)

This class would take a collection of contours (all within the same plane) and generate an appropriate set of polygons. Around half of the code in vtkClipClosedSurface could be put into this class with very little modification. This class would require the contours to be properly connected, so people would be advised to pass their line segments through vtkTriangleFilter before applying it.

This class would be used to do similar things as a constrained vtkDelaunay2D, but would not optimize the triangulation and would not be as sensitive to numerical instability.

vtkCutClosedSurface

This class could be responsible for generating just the cap polygons. On the other hand, this could be provided as an option of the current class.

Guaranteeing a watertight output

Problem

The vtkFeatureEdges filter can check a polydata to see if it is watertight: a watertight surface will have no boundary edges and no non-manifold edges. The vtkClipClosedSurface filter should guarantee in the strictest sense that it will produce a watertight surface if given a watertight surface as an input. Unfortunately, VTK's low-level clipping and triangulation code makes this difficult to achieve.

  1. VTK clipping, as implemented in the vtkCell subclasses, requires the use of a point locator (usually vtkMergePoints) that merges coincident points. Unfortunately, point merging can change connectivity and therefore change topology. Specifically, a manifold data set can, after point merging, become non-manifold.
  2. VTK triangulation, as implemented in vtkPolygon, will leave holes rather than create degenerate or nearly-degenerate triangles. So if the clipping produces polygons that (a) have any 180 degree vertices or (b) have vertices that are very close to each other, then the output will have free edges.

Solution

Unfortunately, the only solution is to do both of these:

  1. Write new clipping code that is specific to this filter, instead of using vtkCell clipping (DONE)
  2. Write new triangulation code, instead of using vtkPolygon::Triangulate()

Another reason not to use vtkPolygon::Triangulate() is performance: vtkPolygon's ear-cut algorithm has a worst-case performance of O(n^2), it was never meant for use on large polygons. A simple locator-based ear-cut algorithm would be O(n) at best and O(n*log n) at worst.

A partial solution is to use a tolerance for creating new points along clipped edges. When an edge is interpolated during the clipping process, if the new point is within tolerance of either of the edge end points, then that endpoint will be used instead of the new point (DONE). Also, if the clipped polygon is found to be degenerate, then the entire polygon can be removed at this stage so that it does not cause problems later.

Reducing the number of triangles produced

Problem

This class uses the clipping/contouring methods in the vtkCell subclasses, which triangulate polys while clipping or contouring them. This produces many more points and cells than are actually needed to describe the geometry. Also, when multiple cut planes are used, these cells can be clipped multiple times, causing a proliferation of small, sharp triangles.

Solution

There are two parts to the solution:

  1. Write new clipping/contouring code that doesn't triangulate the polys while clipping and contouring them. (DONE)
  2. Improve the triangulation of the cap polygons, or even go as far as not triangulating them at all.

The ideal solution is to write clipping/contouring code that can deal with concave polygons. Then, in RequestData where the filter applies each clipping plane in sequence, the code can deal directly with the concave polygons produced by the previous clip-plane, and triangulation of the polygons can be deferred to the very end of RequestData. The result will be that the final output will have the minimum possible number of triangles.

There is another good reason for writing specialized clipping/contouring code: it allows the orientation of the contours to be correctly set from the orientation of the clipped polygons. This will eliminate the need to do an orientation check while creating the cap polygons, and will further guarantee that correct cap polygons are created. For instance, right now vtkClipClosedSurface will create incorrectly-oriented cap polygons if its input is an inside-out surface. But if the orientation of the cap polygons is set from the orientation of the original polygons (instead of being set according to the clip-plane normal like they are now), then cap polygons will be correctly generated even for inside-out surfaces. (DONE)

Handling non-closed inputs

Problem

The algorithm for creating new polygons assumes that all contours are closed. Currently, when a non-closed contour is found, it's two ends are connected in order to close it. This produces some odd results when the filter is given non-closed inputs.

Solution

There are three possible solutions:

  1. Generate an error when non-closed contours are found, but otherwise continue with the current behaviour.
  2. Throw out these non-closed contours. Probably not a good idea.
  3. Add heuristics for optimally joining non-closed contours. (DONE)

The following heuristics have been added for joining the contours:

  1. When joining contours, the new edge must be on the hull of the total point set of all contours.
  2. The new edge must be of the correct sense, i.e. counterclockwise around the hull.
  3. If there are multiple candidate new edges for a loose end, the shortest is chosen.
  4. If there are no new edge candidates for a loose end, then that contour is removed.

Only joins along the hull are permitted, since it is assumed that the loose ends are the result of contours that intersect the bounds of the data set.