[vtkusers] End locations of vtkStreamTracer
Rocco Gasteiger
post at rocco-gasteiger.de
Fri Sep 30 03:39:52 EDT 2011
Hello Tijmen,
I do not know if you question is still relevant but I will give one
suggestion how you can find the end points of your streamlines. The output
of vtkStreamTracer is a vtkPolyDataSet which contains the streamlines as
multiple line cells. You can iterate over each line to retrieve the vertices
of that line. These vertices are ordered beginning with the first vertex in
that line.
However, you must be careful! If you propagate in both directions within the
vtkStreamTracer you have up to two lines for a "whole" streamline. One for
the forward and one for the backward propagation. These two line segments
are not connected! I wrote a function which starts with the first vertex of
the current line. I checked this vertex for two cases. First case: this
vertex is the start point of a streamline with only one segment, that means
only forward or backward propagation was done. Second case: like the first
case but with another vertex at the same position but it is the start point
of the second line segment. I checked this by searching in the neighborhood
of that vertex, with radius 0.0, if another vertex is there. For this
searching I used vtkPointLocator. I connect the two line segments later in
the right order and discard one of the two vertices. Below I give you the
C++ code for that function. I hope that the comments are sufficient. If not
do not hesitate to ask me.
Best regards, Rocco
void UMDStreamlineStyleMapper::sortFlowVertices()
{
inputStreamlineData->BuildCells();
inputStreamlineData->BuildLinks();
vtkSmartPointer<vtkPolyData> sortedStreamlineData =
vtkSmartPointer<vtkPolyData>::New();
vtkSmartPointer<vtkPoints> sortedStreamlinePoints =
vtkSmartPointer<vtkPoints>::New();
vtkSmartPointer<vtkCellArray> sortedStreamlines =
vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkDoubleArray> streamlineScalars =
vtkSmartPointer<vtkDoubleArray>::New();
vtkSmartPointer<vtkPointLocator> pointLocator =
vtkSmartPointer<vtkPointLocator>::New();
pointLocator->SetDataSet(inputStreamlineData);
pointLocator->BuildLocator();
vtkSmartPointer<vtkIdList> composedLinePointIDs =
vtkSmartPointer<vtkIdList>::New();
vtkSmartPointer<vtkIdList> rearrangedLinePointIDs =
vtkSmartPointer<vtkIdList>::New();
vtkSmartPointer<vtkIdList> identicalSeedPointIdsList =
vtkSmartPointer<vtkIdList>::New();
vtkSmartPointer<vtkIdList> currentLinePointIds =
vtkSmartPointer<vtkIdList>::New();
vtkSmartPointer<vtkIdList> incidentStreamlineIdsList =
vtkSmartPointer<vtkIdList>::New();
vtkSmartPointer<vtkIdList> processedSeedpointIds =
vtkSmartPointer<vtkIdList>::New();
progressBar->setFloatValue(0.0);
bool hasScalars = false;
float maxLength = 0.0;
int newPointId = -1;
int numProcessedStreamLines = -1;
if (inputStreamlineData->GetPointData()->GetScalars())
{
streamlineScalars->SetName(inputStreamlineData->GetPointData()->GetScalars()
->GetName());
hasScalars = true;
}
vtkCellArray *unprocessedStreamlines =
inputStreamlineData->GetLines();
if (unprocessedStreamlines == NULL)
{
cerr << "No streamlines available ==> no sorting is performed!"
<< endl;
} else
{
unprocessedStreamlines->InitTraversal();
vtkIdType npts, *pts;
while(unprocessedStreamlines->GetNextCell(npts, pts) && npts >
0)
{
identicalSeedPointIdsList->Reset();
incidentStreamlineIdsList->Reset();
currentLinePointIds->Reset();
rearrangedLinePointIDs->Reset();
composedLinePointIDs->Reset();
numProcessedStreamLines++;
// The first point of a streamline corresponds to its seed
point. We store its id to skip processing a second streamline
// propagation because we already have processed this line
(see below)
int id = pts[0];
if (processedSeedpointIds->IsId(id) != -1)
{
// Continue with the next streamline starting from
a new seed point
continue;
}
// An unprocessed seed point is found. Test if this point
exists twice at the same position.
// The second point is used as seed point for another
streamline (forward or backward propagation).
processedSeedpointIds->InsertNextId(id);
// Get the position of the current seed point
double seedPoint[3];
inputStreamlineData->GetPoint(id, seedPoint);
pointLocator->FindPointsWithinRadius(0.0, seedPoint,
identicalSeedPointIdsList);
int numIdenticalSeedPoints =
identicalSeedPointIdsList->GetNumberOfIds();
int secondId = 0;
int firstId = 0;
switch (numIdenticalSeedPoints)
{
case 1:
// Get all the point ids of the single streamline
(in "identicalSeedPointIdsList" and "incidentStreamlineIdsList"
// are only one entry).
firstId = identicalSeedPointIdsList->GetId(0);
inputStreamlineData->GetPointCells(firstId,
incidentStreamlineIdsList);
inputStreamlineData->GetCellPoints(incidentStreamlineIdsList->GetId(0),
currentLinePointIds);
// Reverse streamline vertex ids because they start
from the seed point.
for (int pointId = 0; pointId <
currentLinePointIds->GetNumberOfIds(); pointId++)
{
// Reverse the ids because they start at the
seed point.
int currentPointId =
currentLinePointIds->GetId(currentLinePointIds->GetNumberOfIds()-pointId-1);
rearrangedLinePointIDs->InsertNextId(currentPointId);
}
break;
case 2:
// Get all the point ids of the first streamline (in
"identicalSeedPointIdsList" are two entries
// and in "incidentStreamlineIdsList" is still one
entry).
inputStreamlineData->GetPointCells(identicalSeedPointIdsList->GetId(0),
incidentStreamlineIdsList);
inputStreamlineData->GetCellPoints(incidentStreamlineIdsList->GetId(0),
currentLinePointIds);
// Reverse streamline vertex ids because they start
from the seed point.
for (int pointId = 0; pointId <
currentLinePointIds->GetNumberOfIds(); pointId++)
{
// Reverse the ids because they start at the
seed point.
int currentPointId =
currentLinePointIds->GetId(currentLinePointIds->GetNumberOfIds()-pointId-1);
rearrangedLinePointIDs->InsertNextId(currentPointId);
}
// Get all the point ids of the second streamline
currentLinePointIds->Reset();
incidentStreamlineIdsList->Reset();
inputStreamlineData->GetPointCells(identicalSeedPointIdsList->GetId(1),
incidentStreamlineIdsList);
secondId = incidentStreamlineIdsList->GetId(0);
inputStreamlineData->GetCellPoints(incidentStreamlineIdsList->GetId(0),
currentLinePointIds);
// We do not reverse the vertex order because we
continue the streamline direction further.
// We also skip the first point of the second
streamline because this point is already inserted
// in the rearranged list.
for (int pointId = 1; pointId <
currentLinePointIds->GetNumberOfIds(); pointId++)
{
// Reverse the ids because they start at the
seed point.
int currentPointId =
currentLinePointIds->GetId(pointId);
rearrangedLinePointIDs->InsertNextId(currentPointId);
}
break;
default:
cerr <<
"UMDStreamlineStyleMapper::sortFlowVertices()==> Invalid number of
streamlines at current seed point!" << endl;
break;
}
float currentLength = 0.0;
for (int pointId = 0; pointId <
rearrangedLinePointIDs->GetNumberOfIds(); pointId++)
{
newPointId++;
int currentPointId =
rearrangedLinePointIDs->GetId(pointId);
if (hasScalars)
{
double scalarValue =
inputStreamlineData->GetPointData()->GetScalars()->GetTuple1(currentPointId)
;
streamlineScalars->InsertNextTuple1(scalarValue);
}
composedLinePointIDs->InsertNextId(newPointId);
sortedStreamlinePoints->InsertNextPoint(inputStreamlineData->GetPoint(curren
tPointId));
if (newPointId > 0)
{
double p0[3];
double p1[3];
sortedStreamlinePoints->GetPoint(newPointId-1,
p0);
sortedStreamlinePoints->GetPoint(newPointId,
p1);
currentLength +=
sqrt(vtkMath::Distance2BetweenPoints(p0, p1));
}
}
sortedStreamlines->InsertNextCell(composedLinePointIDs->GetNumberOfIds(),
composedLinePointIDs->GetPointer(0));
if (currentLength > maxLength)
{
maxLength = currentLength;
}
setMaximumStreamlineLength(maxLength);
// The first half of the progress state consists of
sorting the flow vertices.
// The second half will be done during the generating of
the edge mesh data structure.
progressBar->setFloatValue((float)numProcessedStreamLines/inputStreamlineDat
a->GetNumberOfLines());
} // end while-loop
cout << "sortedPoints->GetNumberOfPoints(): " <<
sortedStreamlinePoints->GetNumberOfPoints() << endl;
//cout << "unsortedPolyLines->GetNumberOfPoints(): " <<
unsortedPolyLines->GetNumberOfPoints() << endl;
//cout << "polyLines->GetNumberOfCells(): " <<
polyLines->GetNumberOfCells() << endl;
//cout << "maxLength: " << maxLength << endl;
sortedStreamlineData->SetPoints(sortedStreamlinePoints);
sortedStreamlineData->SetLines(sortedStreamlines);
sortedStreamlineData->GetPointData()->SetScalars(streamlineScalars);
sortedStreamlineData->Update();
inputStreamlineData = sortedStreamlineData;
vtkOutputSortedStreamlines.setNewOutputBaseFieldObject(sortedStreamlineData)
;
}
}
Von: vtkusers-bounces at vtk.org [mailto:vtkusers-bounces at vtk.org] Im Auftrag
von Tijmen Klein
Gesendet: Dienstag, 27. September 2011 15:10
An: vtkusers at vtk.org
Betreff: Re: [vtkusers] End locations of vtkStreamTracer
A little update on my own question. I am now able to find the end location
of a streamline, as long as I use 1 input point (setStartPosition() instead
of setSource() ). I can find the end point using this trick:
double* stopPosition;
long numPoints = streamer->GetOutput()->GetNumberOfPoints();
stopPosition = streamer->GetOutput()->GetPoint(numPoints-1);
I can then use this stopPosition for the setStartPosition() method.
Unfortunately, this does not work with multiple inputs. Is there anyone who
could help me on that? Ofcouse, I could create an array of vtkStreamTracers,
each one only tracing a single input point. This should work, but feels
really ugly.
Cheers,
Tijmen
On Wed, Sep 21, 2011 at 4:47 PM, Tijmen Klein <T.R.Klein at student.rug.nl>
wrote:
Hi everyone,
I'm trying to get the end locations of the seeds inserted into a
vtkStreamTracer. I insert some particles using a vtkDataSet, and I want to
know the locations of where these particles ended up. How is this possible
with a vtkStreamTracer?
Some context: I'm working with time dependent data, and I want to use the
end locations as the new seed locations for the next timestep. I only want
to visualize the path (using streamtubes) of the current timestep.
Cheers,
Tijmen
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.vtk.org/pipermail/vtkusers/attachments/20110930/fcef910d/attachment.htm>
More information about the vtkusers
mailing list