Commit 251d5129 authored by Berk Geveci's avatar Berk Geveci
Browse files

1. Simplified code 2. Fixed bug: triangle strips with even number of triangles...

1. Simplified code 2. Fixed bug: triangle strips with even number of triangles were not being flipped properly 3. Added support for reflecting around an arbitrary center
parent a997e3b5
......@@ -24,13 +24,14 @@
#include "vtkPointData.h"
#include "vtkUnstructuredGrid.h"
vtkCxxRevisionMacro(vtkReflectionFilter, "1.9");
vtkCxxRevisionMacro(vtkReflectionFilter, "1.10");
vtkStandardNewMacro(vtkReflectionFilter);
//---------------------------------------------------------------------------
vtkReflectionFilter::vtkReflectionFilter()
{
this->Plane = VTK_USE_X_MIN;
this->Plane = USE_X_MIN;
this->Center = 0.0;
}
//---------------------------------------------------------------------------
......@@ -38,6 +39,15 @@ vtkReflectionFilter::~vtkReflectionFilter()
{
}
//---------------------------------------------------------------------------
void vtkReflectionFilter::FlipVector(float tuple[3], int mirrorDir[3])
{
for(int j=0; j<3; j++)
{
tuple[j] *= mirrorDir[j];
}
}
//---------------------------------------------------------------------------
void vtkReflectionFilter::Execute()
{
......@@ -54,8 +64,9 @@ void vtkReflectionFilter::Execute()
float tuple[3];
vtkPoints *outPoints;
float point[3];
float constant;
int ptId, cellId;
float constant[3] = {0.0, 0.0, 0.0};
int mirrorDir[3] = { 1, 1, 1};
int ptId, cellId, j;
vtkGenericCell *cell = vtkGenericCell::New();
vtkIdList *ptIds = vtkIdList::New();
......@@ -91,141 +102,68 @@ void vtkReflectionFilter::Execute()
// Copy reflected points.
switch (this->Plane)
{
case VTK_USE_X_MIN:
constant = 2*bounds[0];
for (i = 0; i < numPts; i++)
{
input->GetPoint(i, point);
ptId =
outPoints->InsertNextPoint(-point[0] + constant, point[1], point[2]);
outPD->CopyData(inPD, i, ptId);
if (inPtVectors)
{
inPtVectors->GetTuple(i, tuple);
tuple[0] = -tuple[0];
outPtVectors->SetTuple(ptId, tuple);
}
if (inPtNormals)
{
inPtNormals->GetTuple(i, tuple);
tuple[0] = -tuple[0];
outPtNormals->SetTuple(ptId, tuple);
}
}
case USE_X_MIN:
constant[0] = 2*bounds[0];
mirrorDir[0] = -1;
break;
case VTK_USE_X_MAX:
constant = 2*bounds[1];
for (i = 0; i < numPts; i++)
{
input->GetPoint(i, point);
ptId =
outPoints->InsertNextPoint(-point[0] + constant, point[1], point[2]);
outPD->CopyData(inPD, i, ptId);
if (inPtVectors)
{
inPtVectors->GetTuple(i, tuple);
tuple[0] = -tuple[0];
outPtVectors->SetTuple(ptId, tuple);
}
if (inPtNormals)
{
inPtNormals->GetTuple(i, tuple);
tuple[0] = -tuple[0];
outPtNormals->SetTuple(ptId, tuple);
}
}
case USE_X_MAX:
constant[0] = 2*bounds[1];
mirrorDir[0] = -1;
break;
case VTK_USE_Y_MIN:
constant = 2*bounds[2];
for (i = 0; i < numPts; i++)
{
input->GetPoint(i, point);
ptId =
outPoints->InsertNextPoint(point[0], -point[1] + constant, point[2]);
outPD->CopyData(inPD, i, ptId);
if (inPtVectors)
{
inPtVectors->GetTuple(i, tuple);
tuple[1] = -tuple[1];
outPtVectors->SetTuple(ptId, tuple);
}
if (inPtNormals)
{
inPtNormals->GetTuple(i, tuple);
tuple[1] = -tuple[1];
outPtNormals->SetTuple(ptId, tuple);
}
}
case USE_X:
constant[0] = this->Center;
mirrorDir[0] = -1;
break;
case VTK_USE_Y_MAX:
constant = 2*bounds[3];
for (i = 0; i < numPts; i++)
{
input->GetPoint(i, point);
ptId =
outPoints->InsertNextPoint(point[0], -point[1] + constant, point[2]);
outPD->CopyData(inPD, i, ptId);
if (inPtVectors)
{
inPtVectors->GetTuple(i, tuple);
tuple[1] = -tuple[1];
outPtVectors->SetTuple(ptId, tuple);
}
if (inPtNormals)
{
inPtNormals->GetTuple(i, tuple);
tuple[1] = -tuple[1];
outPtNormals->SetTuple(ptId, tuple);
}
}
case USE_Y_MIN:
constant[1] = 2*bounds[2];
mirrorDir[1] = -1;
break;
case VTK_USE_Z_MIN:
constant = 2*bounds[4];
for (i = 0; i < numPts; i++)
{
input->GetPoint(i, point);
ptId =
outPoints->InsertNextPoint(point[0], point[1], -point[2] + constant);
outPD->CopyData(inPD, i, ptId);
if (inPtVectors)
{
inPtVectors->GetTuple(i, tuple);
tuple[2] = -tuple[2];
outPtVectors->SetTuple(ptId, tuple);
}
if (inPtNormals)
{
inPtNormals->GetTuple(i, tuple);
tuple[2] = -tuple[2];
outPtNormals->SetTuple(ptId, tuple);
}
}
case USE_Y_MAX:
constant[1] = 2*bounds[3];
mirrorDir[1] = -1;
break;
case VTK_USE_Z_MAX:
constant = 2*bounds[5];
for (i = 0; i < numPts; i++)
{
input->GetPoint(i, point);
ptId =
outPoints->InsertNextPoint(point[0], point[1], -point[2] + constant);
outPD->CopyData(inPD, i, ptId);
if (inPtVectors)
{
inPtVectors->GetTuple(i, tuple);
tuple[2] = -tuple[2];
outPtVectors->SetTuple(ptId, tuple);
}
if (inPtNormals)
{
inPtNormals->GetTuple(i, tuple);
tuple[2] = -tuple[2];
outPtNormals->SetTuple(ptId, tuple);
}
}
case USE_Y:
constant[1] = this->Center;
mirrorDir[1] = -1;
break;
case USE_Z_MIN:
constant[2] = 2*bounds[4];
mirrorDir[2] = -1;
break;
case USE_Z_MAX:
constant[2] = 2*bounds[5];
mirrorDir[2] = -1;
break;
case USE_Z:
constant[2] = this->Center;
mirrorDir[2] = -1;
break;
}
for (i = 0; i < numPts; i++)
{
input->GetPoint(i, point);
ptId =
outPoints->InsertNextPoint( mirrorDir[0]*point[0] + constant[0],
mirrorDir[1]*point[1] + constant[1],
mirrorDir[2]*point[2] + constant[2] );
outPD->CopyData(inPD, i, ptId);
if (inPtVectors)
{
inPtVectors->GetTuple(i, tuple);
this->FlipVector(tuple, mirrorDir);
outPtVectors->SetTuple(ptId, tuple);
}
if (inPtNormals)
{
inPtNormals->GetTuple(i, tuple);
this->FlipVector(tuple, mirrorDir);
outPtNormals->SetTuple(ptId, tuple);
}
}
int numCellPts, j, cellType;
int numCellPts, cellType;
vtkIdType *newCellPts;
vtkIdList *cellPts;
......@@ -243,55 +181,46 @@ void vtkReflectionFilter::Execute()
input->GetCell(i, cell);
numCellPts = cell->GetNumberOfPoints();
cellType = cell->GetCellType();
newCellPts = new vtkIdType[numCellPts];
cellPts = cell->GetPointIds();
for (j = numCellPts-1; j >= 0; j--)
// Triangle strips with even number of triangles have
// to be handled specially. A degenerate triangle is
// introduce to flip all the triangles properly.
if (cellType == VTK_TRIANGLE_STRIP && numCellPts % 2 == 0)
{
numCellPts++;
newCellPts = new vtkIdType[numCellPts];
newCellPts[0] = cellPts->GetId(0) + numPts;
newCellPts[1] = cellPts->GetId(2) + numPts;
newCellPts[2] = cellPts->GetId(1) + numPts;
newCellPts[3] = cellPts->GetId(2) + numPts;
for (j = 4; j < numCellPts; j++)
{
newCellPts[j] = cellPts->GetId(j-1) + numPts;
}
}
else
{
newCellPts[numCellPts-1-j] = cellPts->GetId(j) + numPts;
newCellPts = new vtkIdType[numCellPts];
for (j = numCellPts-1; j >= 0; j--)
{
newCellPts[numCellPts-1-j] = cellPts->GetId(j) + numPts;
}
}
cellId = output->InsertNextCell(cellType, numCellPts, newCellPts);
delete [] newCellPts;
outCD->CopyData(inCD, i, cellId);
if (inCellVectors)
{
inCellVectors->GetTuple(i, tuple);
switch (this->Plane)
{
case VTK_USE_X_MIN:
case VTK_USE_X_MAX:
tuple[0] = -tuple[0];
break;
case VTK_USE_Y_MIN:
case VTK_USE_Y_MAX:
tuple[1] = -tuple[1];
break;
case VTK_USE_Z_MIN:
case VTK_USE_Z_MAX:
tuple[2] = -tuple[2];
break;
}
this->FlipVector(tuple, mirrorDir);
outCellVectors->SetTuple(cellId, tuple);
}
if (inCellNormals)
{
inCellNormals->GetTuple(i, tuple);
switch (this->Plane)
{
case VTK_USE_X_MIN:
case VTK_USE_X_MAX:
tuple[0] = -tuple[0];
break;
case VTK_USE_Y_MIN:
case VTK_USE_Y_MAX:
tuple[1] = -tuple[1];
break;
case VTK_USE_Z_MIN:
case VTK_USE_Z_MAX:
tuple[2] = -tuple[2];
break;
}
this->FlipVector(tuple, mirrorDir);
outCellNormals->SetTuple(cellId, tuple);
}
delete [] newCellPts;
}
cell->Delete();
......
......@@ -27,13 +27,6 @@
#include "vtkDataSetToUnstructuredGridFilter.h"
#define VTK_USE_X_MIN 0
#define VTK_USE_Y_MIN 1
#define VTK_USE_Z_MIN 2
#define VTK_USE_X_MAX 3
#define VTK_USE_Y_MAX 4
#define VTK_USE_Z_MAX 5
class VTK_GRAPHICS_EXPORT vtkReflectionFilter : public vtkDataSetToUnstructuredGridFilter
{
public:
......@@ -42,15 +35,41 @@ public:
vtkTypeRevisionMacro(vtkReflectionFilter, vtkDataSetToUnstructuredGridFilter);
void PrintSelf(ostream &os, vtkIndent indent);
vtkSetClampMacro(Plane, int, 0, 5);
//BTX
enum ReflectionPlane
{
USE_X_MIN = 0,
USE_Y_MIN = 1,
USE_Z_MIN = 2,
USE_X_MAX = 3,
USE_Y_MAX = 4,
USE_Z_MAX = 5,
USE_X = 6,
USE_Y = 7,
USE_Z = 8
};
//ETX
// Description:
// Set the normal of the plane to use as mirror.
vtkSetClampMacro(Plane, int, 0, 8);
vtkGetMacro(Plane, int);
void SetPlaneToXMin() { this->SetPlane(VTK_USE_X_MIN); };
void SetPlaneToYMin() { this->SetPlane(VTK_USE_Y_MIN); };
void SetPlaneToZMin() { this->SetPlane(VTK_USE_Z_MIN); };
void SetPlaneToXMax() { this->SetPlane(VTK_USE_X_MAX); };
void SetPlaneToYMax() { this->SetPlane(VTK_USE_Y_MAX); };
void SetPlaneToZMax() { this->SetPlane(VTK_USE_Z_MAX); };
void SetPlaneToX() { this->SetPlane(USE_X); };
void SetPlaneToY() { this->SetPlane(USE_Y); };
void SetPlaneToZ() { this->SetPlane(USE_Z); };
void SetPlaneToXMin() { this->SetPlane(USE_X_MIN); };
void SetPlaneToYMin() { this->SetPlane(USE_Y_MIN); };
void SetPlaneToZMin() { this->SetPlane(USE_Z_MIN); };
void SetPlaneToXMax() { this->SetPlane(USE_X_MAX); };
void SetPlaneToYMax() { this->SetPlane(USE_Y_MAX); };
void SetPlaneToZMax() { this->SetPlane(USE_Z_MAX); };
// Description:
// If the reflection plane is set to X, Y or Z, this variable
// is use to set the position of the plane.
vtkSetMacro(Center, float);
vtkGetMacro(Center, float);
protected:
vtkReflectionFilter();
~vtkReflectionFilter();
......@@ -58,6 +77,9 @@ protected:
void Execute();
int Plane;
float Center;
void FlipVector(float tuple[3], int mirrorDir[3]);
private:
vtkReflectionFilter(const vtkReflectionFilter&); // Not implemented
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment