[vtkusers] Topology problem in converting polydata to image
donyyo
donyyo at gmail.com
Mon Aug 26 08:20:06 EDT 2013
Hi,
I've deal with this problem by myself. I would like to share my method.
The topology problem posted in my question result from the contour polydata
contains some lines which are not in a circal.
To fill the contour cut by a plane from a polydata, we should know whether a
point is inside or outside the contour. VTK creates a line which is parallel
to axis and end up to the point. If the line intersect with the contour for
an odd number of times, the point is defined to be inside the contour. If it
is even times, the point is in the outside.
I create a simple contour contains a rectangle and four curves (fig. 1 (a)).
The segment image is wrong (fig. 1 (b)).
<http://vtk.1045678.n5.nabble.com/file/n5723025/fig1.png>
A simple method to deal this problem is to remove the lines which don't form
a circal (fig. 2 (a)). Then the image is correct (fig. 2 (b)).
<http://vtk.1045678.n5.nabble.com/file/n5723025/fig2.png>
Hope to useful.
Donny
my code (VS 2005 + VTK 5.8.0) and sample files:
#include "vtkSmartPointer.h"
#include "vtkLineSource.h"
#include "vtkAppendPolyData.h"
#include "vtkCleanPolyData.h"
#include "vtkPolyData.h"
#include "vtkPoints.h"
#include "vtkImageData.h"
#include "vtkPointData.h"
#include "vtkPolyDataToImageStencil.h"
#include "vtkImageStencil.h"
#include "vtkPNGWriter.h"
#include "vtkPolyDataMapper.h"
#include "vtkActor.h"
#include "vtkProperty.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include <iomanip>
#define PRECISION 2
#include <string>
#include <iostream>
using namespace std;
int main(int argc, char *argv[])
{
cout << "---- Initialize point ----" << endl;
vtkSmartPointer<vtkRenderer> ren = vtkSmartPointer<vtkRenderer>::New();
ren->SetBackground(0.1, 0.2, 0.4);
vtkSmartPointer<vtkRenderWindow> renWin =
vtkSmartPointer<vtkRenderWindow>::New();
renWin->AddRenderer(ren);
renWin->SetSize(500, 500);
vtkSmartPointer<vtkRenderWindowInteractor> iren =
vtkSmartPointer<vtkRenderWindowInteractor>::New();
iren->SetRenderWindow(renWin);
string pointFilePath = "F:\\test\\";
string pointFileName = "point2.txt";
string pointFilePathName = pointFilePath+pointFileName;
cout << "-> " << pointFilePathName.c_str() << endl;
ifstream pointFile;
pointFile.open(pointFilePathName.c_str(), ios::in);
if (!pointFile)
{
cerr << "ERROR! Can't open " << pointFilePathName.c_str() << endl;
exit(1);
}
vtkIdType sumPoint = 18;
double **point = new double*[sumPoint];
for (vtkIdType ii = 0; ii < sumPoint; ii++)
point[ii] = new double[3]();
for (vtkIdType ii = 0; ii < sumPoint; ii++)
{
pointFile >> point[ii][0];
pointFile >> point[ii][1];
pointFile >> point[ii][2];
cout << fixed << setprecision(PRECISION) << "point[" << ii << "] =
(" << point[ii][0] << ", " << point[ii][1] << ", " << point[ii][2] << ")" <<
endl;
}
pointFile.close();
string connectionFilePath = pointFilePath;
string connectionFileName = "connection2.txt";
string connectionFilePathName = connectionFilePath+connectionFileName;
cout << "-> " << connectionFilePathName.c_str() << endl;
ifstream connectionFile;
connectionFile.open(connectionFilePathName.c_str(), ios::in);
if (!connectionFile)
{
cerr << "ERROR! Can't open " << connectionFilePathName.c_str() <<
endl;
exit(1);
}
vtkIdType sumLine = 14;
vtkIdType **connection = new vtkIdType*[sumLine];
for (vtkIdType jj = 0; jj < sumLine; jj++)
connection[jj] = new vtkIdType[2]();
for (vtkIdType jj = 0; jj < sumLine; jj++)
{
connectionFile >> connection[jj][0];
connectionFile >> connection[jj][1];
cout << "connection[" << jj << "] = (" << connection[jj][0] << ", "
<< connection[jj][1] << ")" << endl;
}
connectionFile.close();
cout << "---- Get contour ----" << endl;
vtkSmartPointer<vtkAppendPolyData> appendContourPolyData =
vtkSmartPointer<vtkAppendPolyData>::New();
// Construct rectangle.
for (vtkIdType jj = 0; jj < sumLine; jj++)
{
vtkSmartPointer<vtkLineSource> lineSource =
vtkSmartPointer<vtkLineSource>::New();
lineSource->SetPoint1(point[connection[jj][0]]);
lineSource->SetPoint2(point[connection[jj][1]]);
lineSource->Update();
vtkSmartPointer<vtkPolyData> linePolyData =
vtkSmartPointer<vtkPolyData>::New();
linePolyData = lineSource->GetOutput();
appendContourPolyData->AddInput(linePolyData);
appendContourPolyData->Update();
}
vtkSmartPointer<vtkCleanPolyData> cleanContourPolyData =
vtkSmartPointer<vtkCleanPolyData>::New();
cleanContourPolyData->SetInput(appendContourPolyData->GetOutput());
cleanContourPolyData->Update();
cout << "---- Get contourPolyData ----" << endl;
vtkSmartPointer<vtkPolyData> contourPolyData =
vtkSmartPointer<vtkPolyData>::New();
contourPolyData = cleanContourPolyData->GetOutput();
vtkSmartPointer<vtkPolyDataMapper> contourMapper =
vtkSmartPointer<vtkPolyDataMapper>::New();
contourMapper->SetInput(contourPolyData);
vtkSmartPointer<vtkActor> contourActor =
vtkSmartPointer<vtkActor>::New();
contourActor->SetMapper(contourMapper);
contourActor->GetProperty()->SetRepresentationToWireframe();
contourActor->GetProperty()->SetColor(0.0, 1.0, 0.0);
cout << "---- createNewContourPolyData ----" << endl;
vtkSmartPointer<vtkPolyData> newContourPolyData =
vtkSmartPointer<vtkPolyData>::New();
vtkSmartPointer<vtkIdList> cellPointIds =
vtkSmartPointer<vtkIdList>::New();
bool hasUsedOncePoint = true;
while (hasUsedOncePoint == true) // Remove lines which are not in a
circle.
{
int *pointUsedMarker = new
int[contourPolyData->GetNumberOfPoints()]();
for (vtkIdType ii = 0; ii < contourPolyData->GetNumberOfCells();
ii++)
{
contourPolyData->GetCellPoints(ii, cellPointIds);
pointUsedMarker[cellPointIds->GetId(0)]++;
pointUsedMarker[cellPointIds->GetId(1)]++;
}
vtkSmartPointer<vtkAppendPolyData> appendNewContourPolyData =
vtkSmartPointer<vtkAppendPolyData>::New();
vtkSmartPointer<vtkPoints> contourPoints =
vtkSmartPointer<vtkPoints>::New();
contourPoints = contourPolyData->GetPoints();
double cellPoint[2][3];
for (vtkIdType ii = 0; ii < contourPolyData->GetNumberOfCells();
ii++)
{
contourPolyData->GetCellPoints(ii, cellPointIds);
for (vtkIdType jj = 0; jj < 2; jj++)
contourPolyData->GetPoint(cellPointIds->GetId(jj),
cellPoint[jj]);
// If cells contain the points which are used only once, the
cells will not be appended to the new contour polydata.
if (pointUsedMarker[cellPointIds->GetId(0)] >= 2 &&
pointUsedMarker[cellPointIds->GetId(1)] >= 2)
{
vtkSmartPointer<vtkLineSource> lineSource =
vtkSmartPointer<vtkLineSource>::New();
lineSource->SetPoint1(cellPoint[0]);
lineSource->SetPoint2(cellPoint[1]);
lineSource->Update();
appendNewContourPolyData->AddInput(lineSource->GetOutput());
appendNewContourPolyData->Update();
}
}
delete []pointUsedMarker;
vtkSmartPointer<vtkCleanPolyData> cleanNewContourPolyData =
vtkSmartPointer<vtkCleanPolyData>::New();
cleanNewContourPolyData->SetInput(appendNewContourPolyData->GetOutput());
cleanNewContourPolyData->Update();
newContourPolyData->DeepCopy(cleanNewContourPolyData->GetOutput());
int *newPointUsedMarker = new
int[newContourPolyData->GetNumberOfPoints()]();
for (vtkIdType ii = 0; ii < newContourPolyData->GetNumberOfCells();
ii++)
{
newContourPolyData->GetCellPoints(ii, cellPointIds);
newPointUsedMarker[cellPointIds->GetId(0)]++;
newPointUsedMarker[cellPointIds->GetId(1)]++;
}
hasUsedOncePoint = false;
for (vtkIdType jj = 0; jj < newContourPolyData->GetNumberOfPoints();
jj++)
{
if (newPointUsedMarker[jj] < 2)
{
hasUsedOncePoint = true;
break;
}
}
delete []newPointUsedMarker;
contourPolyData->DeepCopy(newContourPolyData);
}
vtkSmartPointer<vtkPolyDataMapper> newContourMapper =
vtkSmartPointer<vtkPolyDataMapper>::New();
newContourMapper->SetInput(newContourPolyData);
vtkSmartPointer<vtkActor> newContourActor =
vtkSmartPointer<vtkActor>::New();
newContourActor->SetMapper(newContourMapper);
newContourActor->GetProperty()->SetRepresentationToWireframe();
newContourActor->GetProperty()->SetColor(1.0, 0.0, 0.0);
cout << "---- Get imageData ----" << endl;
double contourBounds[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
contourPolyData->GetBounds(contourBounds);
cout << "contourBounds = (" << contourBounds[0] << ":" <<
contourBounds[1] << ", " << contourBounds[2] << ":" << contourBounds[3] <<
", " << contourBounds[4] << ":" << contourBounds[5] << ")" << endl;
double imageSpacing[3] = {1.0, 1.0, 1.0};
cout << "imageSpacing = " << imageSpacing[0] << " * " << imageSpacing[1]
<< " * " << imageSpacing[2] << endl;
int imageDimensions[3] = {0, 0, 0};
imageDimensions[0] =
static_cast<int>(contourBounds[1]/imageSpacing[0]+contourBounds[0]/imageSpacing[0]);
imageDimensions[1] =
static_cast<int>(contourBounds[3]/imageSpacing[1]+contourBounds[2]/imageSpacing[1]);
imageDimensions[2] = 1;
cout << "imageDimensions = " << imageDimensions[0] << " * " <<
imageDimensions[1] << " * " << imageDimensions[2] << endl;
int imageExtent[6] = {0, 0, 0, 0, 0, 0};
imageExtent[0] = 0; imageExtent[1] = imageDimensions[0]-1;
imageExtent[2] = 0; imageExtent[3] = imageDimensions[1]-1;
imageExtent[4] = 0; imageExtent[5] = imageDimensions[2]-1;
cout << "imageExtent = (" << imageExtent[0] << ":" << imageExtent[1] <<
", " << imageExtent[2] << ":" << imageExtent[3] << ", " << imageExtent[4] <<
":" << imageExtent[5] << ")" << endl;
double imageOrigin[3] = {0.0, 0.0, 0.0};
imageOrigin[2] = contourBounds[4];
cout << "imageOrigin = (" << imageOrigin[0] << ", " << imageOrigin[1] <<
", " << imageOrigin[2] << ")" << endl;
vtkSmartPointer<vtkImageData> imageData =
vtkSmartPointer<vtkImageData>::New();
imageData->SetSpacing(imageSpacing);
imageData->SetDimensions(imageDimensions);
imageData->SetExtent(imageExtent);
imageData->SetOrigin(imageOrigin);
imageData->SetScalarTypeToUnsignedChar();
imageData->AllocateScalars();
// Fill image with foreground voxels.
for (vtkIdType ii = 0; ii < imageData->GetNumberOfPoints(); ii++)
imageData->GetPointData()->GetScalars()->SetTuple1(ii, 255);
// polygonal data -> image stencil.
vtkSmartPointer<vtkPolyDataToImageStencil> polyDataToImageStencil =
vtkSmartPointer<vtkPolyDataToImageStencil>::New();
polyDataToImageStencil->SetInput(newContourPolyData);
polyDataToImageStencil->SetTolerance(0);
polyDataToImageStencil->SetOutputOrigin(0.0, 0.0, 1.0);
polyDataToImageStencil->SetOutputSpacing(1.0, 1.0, 1.0);
polyDataToImageStencil->SetOutputWholeExtent(imageData->GetExtent());
polyDataToImageStencil->Update();
// Cut the corresponding white image and set the background.
vtkSmartPointer<vtkImageStencil> imageStencil =
vtkSmartPointer<vtkImageStencil>::New();
imageStencil->SetInput(imageData);
imageStencil->SetStencil(polyDataToImageStencil->GetOutput());
imageStencil->ReverseStencilOff();
imageStencil->SetBackgroundValue(0.0);
imageStencil->Update();
cout << "---- Output image ----" << endl;
string outputImageFilePath = "F:\\test\\";
string outputImageFileName = "image.png";
string outputImageFilePathName =
outputImageFilePath+outputImageFileName;
cout << "<- " << outputImageFilePathName.c_str() << endl;
vtkSmartPointer<vtkPNGWriter> pngWriter =
vtkSmartPointer<vtkPNGWriter>::New();
pngWriter->SetInputConnection(imageStencil->GetOutputPort());
pngWriter->SetFileName(outputImageFilePathName.c_str());
pngWriter->Write();
cout << "---- Visualize ----" << endl;
ren->AddActor(newContourActor);
renWin->Render();
iren->Start();
for (vtkIdType ii = 0; ii < sumPoint; ii++)
delete []point[ii];
delete point;
for (vtkIdType jj = 0; jj < sumLine; jj++)
delete []connection[jj];
delete []connection;
return EXIT_SUCCESS;
}
The sample files used in the code contain a point set file (point2.txt) and
connection information file (connection2.txt).
point2.txt
40 40 1
40 80 1
100 80 1
100 40 1
30 40 1
20 60 1
30 80 1
50 90 1
60 100 1
80 100 1
90 90 1
110 80 1
120 60 1
110 40 1
90 30 1
80 20 1
60 20 1
50 30 1
connection2.txt
0 1
1 2
2 3
3 0
4 5
5 6
7 8
8 9
9 10
11 12
12 13
14 15
15 16
16 17
--
View this message in context: http://vtk.1045678.n5.nabble.com/Topology-problem-in-converting-polydata-to-image-tp5722850p5723025.html
Sent from the VTK - Users mailing list archive at Nabble.com.
More information about the vtkusers
mailing list