I wanted to find the area in the black circle on the tooth surface, so the algorithm class vtkSelectPolyData was used in my project.
The idealized result looks like the following picture.

Unfortunately, I got the bad region as this scene.

After reading the source code about vtkSelectPolyData, I found that the points on the circle have to be sorted, or the closed points on the tooth to the black circle will be out-of-order and the selection edges look very messy.

The code about show the selection edges.

    #define vSP vtkSmartPointer
    #define vSPNew(Var, Type) vSP<Type> Var = vSP<Type>::New();

    vtkPolyData *edgeData = selectPd->GetSelectionEdges();
    vSPNew( edgeMapper, vtkPolyDataMapper );
    edgeMapper->SetInputData( edgeData );
    vSPNew( edgeActor, vtkActor );
    edgeActor->SetMapper( edgeMapper );
    edgeActor->GetProperty()->SetColor( 1, 0, 0 );
    renderer->AddActor( edgeActor );

So the points on the black circle must be sorted In clockwise/counterclockwise order.

The interfaces GetCellNeighbors and GetCellEdgeNeighbors were invalid for this polyData because it has just lines.
The first algorithm is to find a neighbor line by current point, then collect another point on the line. Relevant implementation details are in the following code snippet.

    // ===================== create sorted points ============================
    polyData->BuildCells();
    polyData->BuildLinks();
    QVector<vtkIdType> sortedIds;

    vtkIdType lastCellId = 0;
    // the first line
    vtkCell *cell = polyData->GetCell( lastCellId );
    vtkIdList *ids = cell->GetPointIds();
    for( int i = 0; i < ids->GetNumberOfIds(); ++i )
        vtkIdType id = ids->GetId( i );
        sortedIds.push_back( id );
    QVector<bool> cellVisited;
    for( int i = 0; i < polyData->GetNumberOfCells(); ++i )
    {
        cellVisited.push_back( false );
    }

    //polyData->GetCellNeighbors( lastCellId, ids, neighborCellIds );
    //polyData->GetCellEdgeNeighbors( lastCellId, 0, 1, neighborCellIds );

    while ( sortedIds.size() < polyData->GetNumberOfPoints() )
    {
        vtkSPtrNew( neighborCellIds, vtkIdList );
        bool foundCell = false;
        int index = sortedIds.size() - 1;
        while( index >= 0 && !foundCell )
        {
            polyData->GetPointCells( sortedIds[index], neighborCellIds );
            for( int i = 0; i < neighborCellIds->GetNumberOfIds(); ++i )
            {
                vtkIdType cellId = neighborCellIds->GetId( i );
                if( cellVisited[ cellId ] )
                {
                    continue;
                }
                cout << "cellId: " << cellId << endl;
                vtkCell *cell = polyData->GetCell( cellId );
                vtkIdList *ptIds = cell->GetPointIds();
                for( int j = 0; j < ptIds->GetNumberOfIds(); ++j )
                {
                    auto ptId = ptIds->GetId( j );
                    if( !sortedIds.contains( ptId ) )
                    {
                        sortedIds.push_back( ptId );
                        cellVisited[lastCellId] = true;
                        lastCellId = cellId;
                        foundCell = true;
                        break;
                    }
                }
                if( foundCell )
                {
                    break;
                }
            }
            index--;
        }
    }
    // ===================== finished: sorted points ============================

It’s useful for common cases, but there are still bad data sets for special circle.

It seems two polylines construct a circle.

So we need a more powerful algorithm to check both directions for every connection point.

The algorithm class had been showed before, relevant example: Find 3D Model’s Boundary In VTK

Use it in our project:

#include <vtkSmartPointer.h>
#include <vtkProperty.h>
#include <vtkSelectPolyData.h>
#include <vtkSphereSource.h>
#include <vtkClipPolyData.h>
#include <vtkProperty.h>
#include <vtkPolyDataMapper.h>
#include <vtkLODActor.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkRenderWindow.h>
#include <vtkRenderer.h>
#include <vtkDataSet.h>
#include <vtkXMLPolyDataReader.h>
#include <vtkTextSource.h>
#include <QString>
#include <vtkCoordinate.h>
#include <vtkTransformPolyDataFilter.h>
#include <vtkPolyDataMapper2D.h>
#include <vtkActor2D.h>
#include <vtkTransform.h>
#include <vtkNamedColors.h>
#include <vtkCell.h>
#include <QVector>
#include <vtkIdList.h>
#include <vtkSphereSource.h>
#include <vtkFollower.h>
#include <vtkPolyDataConnectivityFilter.h>
#include <vtkVectorText.h>
#include <QMap>
#include <QSet>
#include <vtkCellData.h>

#include "ConnectedEdgeFilter.hpp"

int main(int, char *[])
{
    vSPNew( reader, vtkXMLPolyDataReader );
    reader->SetFileName( "../intersectCircle3.vtp" );
    reader->Update();

    auto *polyData = reader->GetOutput();

    vSPNew( connectivityF, vtkPolyDataConnectivityFilter );
    connectivityF->SetInputData( polyData );
    connectivityF->SetExtractionModeToLargestRegion();
    connectivityF->Update();

    polyData = connectivityF->GetOutput();

    vSPNew( mapper, vtkPolyDataMapper );
    mapper->SetInputData( polyData );

    vSPNew( actor, vtkActor );
    actor->SetMapper( mapper );

    vSPNew( renderer, vtkRenderer );
    renderer->AddActor(actor);
    renderer->SetBackground( 0, 0, 0 );

    vSPNew( renderWindow, vtkRenderWindow );
    renderWindow->AddRenderer( renderer );

    vSPNew( renderWindowInteractor, vtkRenderWindowInteractor );
    renderWindowInteractor->SetRenderWindow( renderWindow );

    // ===================== connectivity filter ============================
    ConnectedEdgeFilter* connectFilter = new ConnectedEdgeFilter;
    connectFilter->Initialise( polyData->GetLines() );
    connectFilter->HandleEdges();

    vSP<vtkPoints> connectPoints = vSP<vtkPoints>::New();
    vSP<vtkIdList> longestList = connectFilter->GetLongestList();
    for (vtkIdType i = 0; i < longestList->GetNumberOfIds(); i++)
    {
        connectPoints->InsertNextPoint( polyData->GetPoint(longestList->GetId(i)) );
    }
    delete connectFilter;
    // ===================== Finish: connectivity filter ============================

    for( int i = 0; i < connectPoints->GetNumberOfPoints(); ++i )
    {
        // text 2D
        vtkSmartPointer<vtkTextSource> text2D =
                vtkSmartPointer<vtkTextSource>::New();
        text2D->SetText( QString::number( i ).toStdString().c_str() );

        vtkSmartPointer<vtkTransform> text2DTransform =
                vtkSmartPointer<vtkTransform>::New();
        double *center = connectPoints->GetPoint( i );
        text2DTransform->Translate( center[0], center[1], center[2] );
        text2DTransform->Scale( 0.003, 0.003, 0.003 );

        vtkSmartPointer<vtkTransformPolyDataFilter> text2DDataFilter =
                vtkSmartPointer<vtkTransformPolyDataFilter>::New();
        text2DDataFilter->SetTransform( text2DTransform );
        text2DDataFilter->SetInputConnection( text2D->GetOutputPort() );

        vtkSmartPointer<vtkCoordinate> coords =
                vtkSmartPointer<vtkCoordinate>::New();
        coords->SetCoordinateSystemToWorld();

        vtkSmartPointer<vtkPolyDataMapper2D> text2DMapper =
                vtkSmartPointer<vtkPolyDataMapper2D>::New();
        text2DMapper->SetInputConnection( text2DDataFilter->GetOutputPort() );
        text2DMapper->SetTransformCoordinate( coords );

        vtkSmartPointer<vtkActor2D> text2DActor =
                vtkSmartPointer<vtkActor2D>::New();
        text2DActor->SetMapper( text2DMapper );

        renderer->AddActor( text2DActor );
    }

    renderer->ResetCamera();
    renderWindow->Render();
    renderWindowInteractor->Start();
    return 0;
}

Categories: VTK

1 Comment

vtkSelectPolyData可能存在bug ? | | theArcticOcean · 2020年6月14日 at pm9:26

[…] VTK – Sort The Points On The Circle VTK – Expand Sample Area From A […]

Leave a Reply

Your email address will not be published. Required fields are marked *

You cannot copy content of this page