[VTK]基於VTK的任意平面切割

// 先貼碼 以后再。。。

// 切割介紹
// 對於一個模型的切割需要怎么辦呢，想想切西瓜就知道，首先需要有一個模型、然后有一個切割平面
// 接着對於每個切割操作來更新模型，這樣就可以得到切割的效果了
#include "vtkPlanes.h"
#include "vtkProperty.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkVolume.h"
#include "vtkVolumeProperty.h"
#include "vtkXMLImageDataReader.h" 

#include "vtkContourFilter.h"
#include "vtkSmartPointer.h"
#include "vtkPolyDataNormals.h"
#include "vtkPolyDataMapper.h"
#include "vtkActor.h"
#include "vtkOutlineFilter.h"
#include "vtkStripper.h"

#include "vtkSmoothPolyDataFilter.h"

#include <vtkSphereSource.h>
#include <vtkImagePlaneWidget.h>
#include <vtkInteractorStyleTrackballActor.h>
#include <vtkInteractorStyleTrackballCamera.h>
#include "vtkActor.h"
#include "vtkImageFlip.h"
#include "vtkImageResample.h"
#include "vtkImageViewer.h"
#include "vtkConeSource.h"
#include "vtkBoxWidget.h"
#include "vtkTransform.h"

#include "VTKReBuild.h"
#include <vtkLineSource.h>
#include <vtkDataSetMapper.h>
#include <vector>
#include <OpenCV243.h>

#include "vtkPolyDataWriter.h"
#include "vtkPolyDataReader.h"

using namespace std;
using namespace cv;

class BuildVTKWidgetCall : public vtkCommand
{
public:

 static BuildVTKWidgetCall *New()
 {
  return new BuildVTKWidgetCall;
 }
public:
 virtual void Execute(vtkObject *caller, unsigned long eventId, void *callData)
 {
     vtkImplicitPlaneWidget *pWidget = reinterpret_cast<vtkImplicitPlaneWidget*>(caller);
     if (pWidget)
     {
         vtkSmartPointer<vtkPlane> planeNew = vtkPlane::New();
         pWidget->GetPlane(planeNew);

         cliper->SetClipFunction(planeNew);
         cliper->Update();

         vtkSmartPointer<vtkPolyData> clipedData = vtkPolyData::New();
         clipedData->DeepCopy(cliper->GetOutput());

         vtkSmartPointer<vtkPolyDataMapper> coneMapper = vtkPolyDataMapper::New();
         coneMapper->SetInput(clipedData);
         coneMapper->ScalarVisibilityOff();    
         actor->SetMapper(coneMapper);
     }
 }
 void setCliper(vtkSmartPointer<vtkClipPolyData> other){cliper = other;}
 void setPlane(vtkSmartPointer<vtkPlane> other){pPlane = other;}
 void setActor(vtkSmartPointer<vtkActor> other){actor = other;}
private:
     vtkSmartPointer<vtkPlane> pPlane;
     vtkSmartPointer<vtkActor> actor;
     vtkSmartPointer<vtkClipPolyData> cliper;
};

void build3DView()
{
    vtkSmartPointer<vtkRenderer> aRenderer =
        vtkSmartPointer<vtkRenderer>::New();
    vtkSmartPointer<vtkRenderWindow> renWin =
        vtkSmartPointer<vtkRenderWindow>::New();
    renWin->AddRenderer(aRenderer);

    vtkSmartPointer<vtkRenderWindowInteractor> iren =
        vtkSmartPointer<vtkRenderWindowInteractor>::New();
    iren->SetRenderWindow(renWin);

    vtkSmartPointer<vtkJPEGReader> dicomReader =
        vtkSmartPointer<vtkJPEGReader>::New();  

    dicomReader->SetFilePrefix("C:/Users/DawnWind/Desktop/000/");
    dicomReader->SetFilePattern("%s%d.jpg");
    dicomReader->SetDataByteOrderToLittleEndian();
    dicomReader->SetDataSpacing(1, 1, 1.4); 
    dicomReader->SetFileNameSliceSpacing(1); 
    dicomReader->SetDataExtent(0, 209, 0, 209, 0, 29);
    dicomReader->Update();  

    vtkSmartPointer<vtkContourFilter> skinExtractor =
        vtkSmartPointer<vtkContourFilter>::New();
    skinExtractor->SetInputConnection(dicomReader->GetOutputPort());
    skinExtractor->SetValue(0, 100);    //值越大，保留的部分越少。
#ifdef TEST
    // 這里有用到OpenCV 如果用戶沒有安裝OpenCV那么請將與之有關的刪除
    // 這里一定要update不然下面的getpoints之類是無法取得數據的
    // 這就是在http://www.cnblogs.com/dawnWind/archive/2013/02/17/3D_06.html 里提到的jpegReader->Update()
    skinExtractor->Update();
    auto data = skinExtractor->GetOutput();
    auto points = data->GetPoints();
    auto pSize = points->GetNumberOfPoints();
    vector<Point3d> pointsGroup;
    Mat newMat = Mat::zeros(210, 210, CV_8UC1);
    int matStep = newMat.step;
    auto matData = newMat.data;
    Point2d center;
    for (int i = 0; i < pSize; i++)
    {
        double point[3];
        points->GetPoint(i, point);
        Point3d p1;
        p1.x = (point[0]);
        p1.y = (point[1]);
        p1.z = (point[2]);
        *(matData + (int)point[0] + (int)point[1] * matStep) = 255;
        pointsGroup.push_back(p1);
        center.x += (int)point[0];
        center.y += (int)point[1];
    }
    center.x /= pSize;
    center.y /= pSize;
    imshow("mat", newMat);
    //Mat dst0;
    //flip(newMat, dst0, 0);
    //imshow("dst0", dst0);
    //Mat dst1;
    //flip(newMat, dst1, 1);
    //imshow("dst1", dst1);
    //Mat dstn1;
    //flip(newMat, dstn1, -1);
    //imshow("dstn1", dstn1);
    //imwrite("a.jpg", newMat);
    // 圖像本身是與原始圖像成某軸對稱因此不能在原圖中找中心點
    waitKey();
#endif
    /**做平滑處理**/
    vtkSmartPointer<vtkSmoothPolyDataFilter> smooth = vtkSmoothPolyDataFilter::New();
    smooth->SetInput( skinExtractor->GetOutput());
    smooth->SetNumberOfIterations( 100 );

    //重新計算法向量
    vtkSmartPointer<vtkPolyDataNormals> skinNormals =
        vtkSmartPointer<vtkPolyDataNormals>::New();

    skinNormals->SetInputConnection(smooth->GetOutputPort());
    skinNormals->SetFeatureAngle(60.0);        

    vtkSmartPointer<vtkStripper> skinStripper =        //create triangle strips and/or poly-lines 為了更快的顯示速度
        vtkSmartPointer<vtkStripper>::New();
    skinStripper->SetInputConnection(skinNormals->GetOutputPort()); 

    vtkSmartPointer<vtkPolyDataMapper> skinMapper =
        vtkSmartPointer<vtkPolyDataMapper>::New();
    skinMapper->SetInputConnection(skinStripper->GetOutputPort());
    skinMapper->ScalarVisibilityOff();    //這樣不會帶顏色

    vtkSmartPointer<vtkActor> skin =
        vtkSmartPointer<vtkActor>::New();
    skin->SetMapper(skinMapper); 

    vtkSmartPointer<vtkOutlineFilter> outlineData =
        vtkSmartPointer<vtkOutlineFilter>::New();
    outlineData->SetInputConnection(dicomReader->GetOutputPort());

    vtkSmartPointer<vtkPolyDataMapper> mapOutline =
        vtkSmartPointer<vtkPolyDataMapper>::New();
    mapOutline->SetInputConnection(outlineData->GetOutputPort());

    vtkSmartPointer<vtkActor> outline =
        vtkSmartPointer<vtkActor>::New();
    outline->SetMapper(mapOutline);
    outline->GetProperty()->SetColor(0, 0, 0);
    
    vtkSmartPointer<vtkCamera> aCamera =
        vtkSmartPointer<vtkCamera>::New();
    aCamera->SetViewUp (0, 0, -1);
    aCamera->SetPosition (0, 1, 0);
    aCamera->SetFocalPoint (0, 0, 0);
    aCamera->ComputeViewPlaneNormal();
    aCamera->Azimuth(30.0);
    aCamera->Elevation(30.0);
    aCamera->Dolly(1.5);  

    aRenderer->AddActor(outline);
    aRenderer->AddActor(skin);
    aRenderer->SetActiveCamera(aCamera);
    aRenderer->ResetCamera ();
    aRenderer->SetBackground(.2, .3, .4);
    aRenderer->ResetCameraClippingRange ();

    renWin->SetSize(640, 480);

    vtkSmartPointer<vtkInteractorStyleTrackballCamera> style = 
    vtkSmartPointer<vtkInteractorStyleTrackballCamera>::New();
    iren->SetInteractorStyle( style );
      /////////設置截面
      vtkSmartPointer<vtkClipPolyData> cliper = vtkClipPolyData::New();
      cliper->SetInput(skinStripper->GetOutput());

      // 此平面box可以通過右鍵來進行放大縮小處理（只有當鼠標控制區域只有切割體才單一有效）
      vtkSmartPointer<vtkImplicitPlaneWidget> implicitPlaneWidget = vtkImplicitPlaneWidget::New();
      implicitPlaneWidget->SetInteractor(iren);
      implicitPlaneWidget->SetPlaceFactor(1.25);
      //initially position the widget 
      implicitPlaneWidget->SetInput(skinStripper->GetOutput());
      implicitPlaneWidget->PlaceWidget();
      
      //////Render2
      vtkSmartPointer<vtkActor> coneSkinActor = vtkActor::New();
      coneSkinActor->SetMapper( skinMapper );
      
      vtkSmartPointer<vtkRenderer> rRenderer = 
      vtkSmartPointer<vtkRenderer>::New();
      rRenderer->SetBackground( 0.2, 0.3, 0.5 );
      rRenderer->SetViewport(0.5, 0.0, 1.0, 1.0);
      rRenderer->AddActor(coneSkinActor);

      vtkSmartPointer<BuildVTKWidgetCall> pCall = BuildVTKWidgetCall::New();
      pCall->setActor(coneSkinActor);
      pCall->setCliper(cliper);

      renWin->AddRenderer(rRenderer);
      ///////

    implicitPlaneWidget->AddObserver(vtkCommand::EndInteractionEvent, pCall);
    implicitPlaneWidget->On();
    // Render
    renWin->Render();
    // Initialize the event loop and then start it.
    iren->Initialize();
    
    iren->Start();
}