Search code examples
c++image-processingvtkmedical

How to crop imagedata after a registration?


This is a question about vtk program.

I'm making a rigid registration between 3D MRI and 3D US image. The size of US data is much bigger than MRI. First, I take some landmarks on the same organ of both images.Then, I use "vtkLandmarkTransform" to get a transform between these landmark points. At last I apply the transform to MRI data by using "vtkImageReslice". So I get a huge MRI data.

How can I crop the huge data to make the new MRI data have the same size as US image and show the same organ position as US image?

Hers's my code

vtkSmartPointer<vtkPoints> sourcePoints = vtkSmartPointer<vtkPoints>::New();//Input source points
for ( int t = 0; t < 6; t++ )
{
    float sourcePoint[3] = {inital_points_mri[t][0], inital_points_mri[t][1], inital_points_mri[t][2]};
    sourcePoints->InsertNextPoint(sourcePoint);
}

vtkSmartPointer<vtkPoints> targetPoints = vtkSmartPointer<vtkPoints>::New();//Input target points
for ( int t = 0; t < 6; t++ )
{
    float targetPoint[3] = {inital_points_us[t][0], inital_points_us[t][1], inital_points_us[t][2]};
    targetPoints->InsertNextPoint(targetPoint);
}

vtkSmartPointer<vtkLandmarkTransform> landmarkTransform = vtkSmartPointer<vtkLandmarkTransform>::New();//get landmark transform
landmarkTransform->SetSourceLandmarks(sourcePoints);
landmarkTransform->SetTargetLandmarks(targetPoints);
landmarkTransform->SetModeToSimilarity ();
landmarkTransform->Update();

vtkSmartPointer<vtkImageReslice> transform2 = vtkSmartPointer<vtkImageReslice>::New(); // Apply transform
transform2->SetInput(MRIimagedata); // Input MRI image
transform2->AutoCropOutputOn();  // I'm not sure I need it or not?
// transform2->SetOutputExtent(0,499,0,489,0,358); // set the extent as US data, but can't get the same organ position
landmarkTransform->Inverse();
transform2->SetResliceTransform(landmarkTransform);
transform2->Update();
vtkImageData* transformImage = transform2->GetOutput(); 

Thank you for your attention.


Solution

  • Sorry, guys. I made a stupid mistake in this part.

    I take landmarks of MRI on a coordinate with origin (0, 0, 0). Then I apply the transform to MRI data but forget that the offset of MRI data is not (0, 0, 0). That why I can't get the right result.

    Back to my question: How to crop image data by using VTK? The answer is very simple: Set output origin and extend parameters by using following code:

    vtkSmartPointer<vtkImageReslice> transform2 = vtkSmartPointer<vtkImageReslice>::New();
    transform2->SetInput(MRIimagedata);
    // transform2->AutoCropOutputOn();  //It's not necessary 
    landmarkTransform1->Inverse();
    transform2->SetResliceTransform(landmarkTransform1);
    transform2->SetOutputOrigin(0,0,0);  //set the origin of new image data. It's where we start our crop
    transform2->SetOutputExtent(0,499,0,489,0,358);  //set the range to crop from big data
    

    If the transform is correct, you will get the right result.

    Here's my code

    // read MRI image data
    vtkSmartPointer<vtkMetaImageReader> reader1 =
        vtkSmartPointer<vtkMetaImageReader>::New();
    reader1->SetFileName("MRI.mhd");
    reader1->Update();
    vtkImageData* MRIimagedata = reader1->GetOutput();
    double* origin_MRI = MRIimagedata->GetOrigin();
    
    // generate landmark transforms
    vtkSmartPointer<vtkPoints> sourcePoints1 = vtkSmartPointer<vtkPoints>::New();//Input landmarks on MRI data, should add the offset of MRI data
    for ( int t = 0; t < 6; t++ )
    {
        float sourcePoint1[3] = {inital_points_mri[t][0]+origin_MRI[0], inital_points_mri[t][1]+origin_MRI[1], inital_points_mri[t][2]+origin_MRI[2]};
        sourcePoints1->InsertNextPoint(sourcePoint1);
    }
    
    vtkSmartPointer<vtkPoints> targetPoints1 = vtkSmartPointer<vtkPoints>::New();//input landmarks on US data, the offset of US is (0,0,0)
    for ( int t = 0; t < 6; t++ )
    {
        float targetPoint1[3] = {inital_points_us[t][0], inital_points_us[t][1], inital_points_us[t][2]};
        targetPoints1->InsertNextPoint(targetPoint1);
    }
    
    vtkSmartPointer<vtkLandmarkTransform> landmarkTransform1 = vtkSmartPointer<vtkLandmarkTransform>::New();
    landmarkTransform1->SetSourceLandmarks(sourcePoints1);
    landmarkTransform1->SetTargetLandmarks(targetPoints1);
    landmarkTransform1->SetModeToSimilarity();
    landmarkTransform1->Update();
    
    vtkSmartPointer<vtkImageReslice> transform2 = vtkSmartPointer<vtkImageReslice>::New();
    transform2->SetInput(MRIimagedata);
    // transform2->AutoCropOutputOn();
    landmarkTransform1->Inverse();
    transform2->SetResliceTransform(landmarkTransform1);
    transform2->SetOutputOrigin(0,0,0);
    transform2->SetOutputExtent(0,499,0,489,0,358);
    transform2->Update();
    vtkImageData* transformImage = transform2->GetOutput();