3D реконструкция поверхности сосуда

Question

У меня есть объемный объем ультразвуковой объемной ультразвуковой объем, содержащий один сосуд, и я пытаюсь восстановить поверхность сосуда. 3D-объем создается из стека 2D-изображений/B-сканов, а контур сосуда в каждом B-сканировании сегментирован; то есть я имею эллипс, представляющий контур сосуда в каждом B-сканировании в объеме. Я попытался восстановить контур судна, следуя примеру VTK «GenerateModelsFromLabels.cxx» (http://www.vtk.org/Wiki/VTK/Examples/Cxx/Medical/GenerateModelsFromLabels). Однако результат - это не гладкая поверхность от одного кадра к другому, как я бы надеялся на это. Он прерывистый и нерегулярный, и поверхность не соединяет контуры сосуда между двумя соседними рамками в объеме, если смещение между эллипсами велико. В моем подходе я в основном использовал DiscreteMarchingCubes -> WindowedSincPolyDataFilter -> GeometryFilter.

Я играл с полосой пропускания, smoothingIterations и featureAngle параметров, и я был в состоянии получить лучший следующий результат:

Как вы можете видеть , это не гладкая непрерывная поверхность с множеством неинтерполированных «дыр» между соседними кадрами, но все в порядке. Может ли это быть лучше? Я также попытался использовать триангуляцию 3D Delaunay, но это только дало мне выпуклый корпус, который не был результатом, который я ожидал. Я хотел бы знать, есть ли лучший подход к восстановлению поверхности, которая следует за контуром сосуда от одного B-сканирования до следующего в объеме?

Минимальный рабочий пример приведен ниже:

vtkSmartPointer<vtkImageData> vesselVolume = 
    vtkSmartPointer<vtkImageData>::New(); 

int totalImages = 210; 

for (int z = 0; z < totalImages; z++) 
{ 
    std::string strFile = "E:/datasets/vasc/rendering/contour/" + std::to_string(z + 1) + ".png"; 

    cv::Mat im = cv::imread(strFile, CV_LOAD_IMAGE_GRAYSCALE); 

    if (z == 0) 
    { 
     vesselVolume->SetExtent(0, im.cols, 0, im.rows, 0, totalImages - 1); 
     vesselVolume->SetSpacing(1, 1, 1); 
     vesselVolume->SetOrigin(0, 0, 0); 
     vesselVolume->AllocateScalars(VTK_UNSIGNED_CHAR, 0); 
    } 

    std::vector<cv::Point2i> locations; // output, locations of non-zero pixels 
    cv::findNonZero(im, locations); 

    for (int nzi = 0; nzi < locations.size(); nzi++) 
    { 
     unsigned char* pixel = static_cast<unsigned char*>(vesselVolume->GetScalarPointer(locations[nzi].x, locations[nzi].y, z)); 
     pixel[0] = 255; 
    } 
} 

vtkSmartPointer<vtkDiscreteMarchingCubes> discreteCubes = 
    vtkSmartPointer<vtkDiscreteMarchingCubes>::New(); 

discreteCubes->SetInputData(vesselVolume); 
discreteCubes->GenerateValues(1, 255, 255); 
discreteCubes->ComputeNormalsOn(); 

vtkSmartPointer<vtkWindowedSincPolyDataFilter> smoother = 
    vtkSmartPointer<vtkWindowedSincPolyDataFilter>::New(); 

unsigned int smoothingIterations = 10; 
double passBand = 2; 
double featureAngle = 360.0; 

smoother->SetInputConnection(discreteCubes->GetOutputPort()); 
smoother->SetNumberOfIterations(smoothingIterations); 
smoother->BoundarySmoothingOff(); 
//smoother->FeatureEdgeSmoothingOff(); 
smoother->FeatureEdgeSmoothingOn(); 
smoother->SetFeatureAngle(featureAngle); 
smoother->SetPassBand(passBand); 
smoother->NonManifoldSmoothingOn(); 
smoother->BoundarySmoothingOn(); 
smoother->NormalizeCoordinatesOn(); 
smoother->Update(); 

vtkSmartPointer<vtkThreshold> selector = 
    vtkSmartPointer<vtkThreshold>::New(); 

selector->SetInputConnection(smoother->GetOutputPort()); 
selector->SetInputArrayToProcess(0, 0, 0, 
    vtkDataObject::FIELD_ASSOCIATION_CELLS, 
    vtkDataSetAttributes::SCALARS); 

vtkSmartPointer<vtkMaskFields> scalarsOff = 
    vtkSmartPointer<vtkMaskFields>::New(); 

// Strip the scalars from the output 
scalarsOff->SetInputConnection(selector->GetOutputPort()); 
scalarsOff->CopyAttributeOff(vtkMaskFields::POINT_DATA, 
    vtkDataSetAttributes::SCALARS); 
scalarsOff->CopyAttributeOff(vtkMaskFields::CELL_DATA, 
    vtkDataSetAttributes::SCALARS); 

vtkSmartPointer<vtkGeometryFilter> geometry = 
    vtkSmartPointer<vtkGeometryFilter>::New(); 

geometry->SetInputConnection(scalarsOff->GetOutputPort()); 
geometry->Update(); 

vtkSmartPointer<vtkPolyDataMapper> mapper = 
    vtkSmartPointer<vtkPolyDataMapper>::New(); 
mapper->SetInputConnection(geometry->GetOutputPort()); 
mapper->ScalarVisibilityOff(); 
mapper->Update(); 

vtkSmartPointer<vtkRenderWindow> renderWindow = 
    vtkSmartPointer<vtkRenderWindow>::New(); 

vtkSmartPointer<vtkRenderWindowInteractor> renderWindowInteractor = 
    vtkSmartPointer<vtkRenderWindowInteractor>::New(); 
renderWindowInteractor->SetRenderWindow(renderWindow); 

    vtkSmartPointer<vtkRenderer> renderer = 
    vtkSmartPointer<vtkRenderer>::New(); 

renderWindow->AddRenderer(renderer); 
renderer->SetBackground(.2, .3, .4); 

vtkSmartPointer<vtkActor> actor = 
    vtkSmartPointer<vtkActor>::New(); 
actor->SetMapper(mapper); 
renderer->AddActor(actor); 

renderer->ResetCamera(); 

renderWindow->Render(); 
renderWindowInteractor->Start();   

Answer 1

Предполагая, что ваша проблема рука дрожала между ломтиками, один из возможных способов улучшить свой результат, чтобы применить ломтик порезать регистрации. Попробуйте использовать ImageJ. Используйте трансляции между срезами, чтобы также преобразовать ваши помеченные изображения. Затем запустите преобразованные изображения меток через ваш текущий конвейер.