From: tbaudier Date: Mon, 14 Dec 2015 13:14:03 +0000 (+0100) Subject: Debug Rigid Registration X-Git-Tag: v1.4.0~83^2 X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=commitdiff_plain;h=c86c03213cb8bb9f5c4830d1f96e792a13a34efe;p=clitk.git Debug Rigid Registration With & Without Overlay/Fusion/VF --- diff --git a/vv/vvSlicer.cxx b/vv/vvSlicer.cxx index 8eab74a..63d4699 100644 --- a/vv/vvSlicer.cxx +++ b/vv/vvSlicer.cxx @@ -1080,23 +1080,18 @@ void vvSlicer::UpdateDisplayExtent() #if VTK_MAJOR_VERSION <= 5 input->UpdateInformation(); #else - int extent[6]; mRegisterExtent = mImageReslice->GetOutputInformation(0)->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT()); - copyExtent(mRegisterExtent, extent); - - // Make sure that the required part image has been computed - extent[SliceOrientation*2] = Slice; - extent[SliceOrientation*2+1] = Slice; - - mImageReslice->SetUpdateExtent(extent); #endif this->SetSlice( this->GetSlice() ); //SR: make sure the update let the slice in extents // Local copy of extent int w_ext[6]; - int* ext = GetExtent(); + int* ext = mImageReslice->GetOutputInformation(0)->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT()); copyExtent(ext, w_ext); - + if (mUseReducedExtent) { + copyExtent(mReducedExtent, w_ext); + } + // Set slice value w_ext[ this->SliceOrientation*2 ] = this->Slice; @@ -1105,29 +1100,6 @@ void vvSlicer::UpdateDisplayExtent() // Image actor this->ImageActor->SetVisibility(mImageVisibility); this->ImageActor->SetDisplayExtent(w_ext); -#if VTK_MAJOR_VERSION >= 6 - vtkSmartPointer mapperOpenGL= vtkSmartPointer::New(); - try { - mapperOpenGL = dynamic_cast(GetImageActor()->GetMapper()); - } catch (const std::bad_cast& e) { - std::cerr << e.what() << std::endl; - std::cerr << "Conversion error" << std::endl; - return; - } - if (mFirstSetSliceOrientation) { - copyExtent(ext, mRegisterExtent); - } else { - int w_croppingRegion[6]; - if (mUseReducedExtent) { - copyExtent(mReducedExtent, w_croppingRegion); - } else { - copyExtent(mRegisterExtent, w_croppingRegion); - } - w_croppingRegion[ this->SliceOrientation*2 ] = this->Slice; - w_croppingRegion[ this->SliceOrientation*2+1 ] = this->Slice; - mapperOpenGL->SetCroppingRegion(w_croppingRegion); - } -#endif #if VTK_MAJOR_VERSION >= 6 || (VTK_MAJOR_VERSION >= 5 && VTK_MINOR_VERSION >= 10) // Fix for bug #1882 @@ -1136,10 +1108,9 @@ void vvSlicer::UpdateDisplayExtent() // Overlay image actor if (mOverlay && mOverlayVisibility) { - //mOverlayMapper->GetOutput()->Print(cout); AdjustResliceToSliceOrientation(mOverlayReslice); int overExtent[6]; - this->ConvertImageToImageDisplayExtent(input, w_ext, mOverlayReslice->GetOutput(), overExtent); + this->ConvertImageToImageDisplayExtent(mImageReslice->GetOutputInformation(0), w_ext, mOverlayReslice->GetOutput(), overExtent); #if VTK_MAJOR_VERSION <= 5 bool out = ClipDisplayedExtent(overExtent, mOverlayMapper->GetInput()->GetWholeExtent()); #else @@ -1159,7 +1130,7 @@ void vvSlicer::UpdateDisplayExtent() if (mFusion && mFusionVisibility) { AdjustResliceToSliceOrientation(mFusionReslice); int fusExtent[6]; - this->ConvertImageToImageDisplayExtent(input, w_ext, mFusionReslice->GetOutput(), fusExtent); + this->ConvertImageToImageDisplayExtent(mImageReslice->GetOutputInformation(0), w_ext, mFusionReslice->GetOutput(), fusExtent); #if VTK_MAJOR_VERSION <= 5 bool out = ClipDisplayedExtent(fusExtent, mFusionMapper->GetInput()->GetWholeExtent()); #else @@ -1199,7 +1170,7 @@ void vvSlicer::UpdateDisplayExtent() #else //this->UpdateInformation(); #endif - this->ConvertImageToImageDisplayExtent(input, w_ext, mVF->GetVTKImages()[0], vfExtent); + this->ConvertImageToImageDisplayExtent(mImageReslice->GetOutputInformation(0), w_ext, mVF->GetVTKImages()[0], vfExtent); #if VTK_MAJOR_VERSION <= 5 bool out = ClipDisplayedExtent(vfExtent, mVOIFilter->GetInput()->GetWholeExtent()); #else @@ -1258,13 +1229,16 @@ void vvSlicer::UpdateDisplayExtent() //---------------------------------------------------------------------------- //---------------------------------------------------------------------------- -void vvSlicer::ConvertImageToImageDisplayExtent(vtkImageData *sourceImage, const int sourceExtent[6], +void vvSlicer::ConvertImageToImageDisplayExtent(vtkInformation *sourceImage, const int sourceExtent[6], vtkImageData *targetImage, int targetExtent[6]) { //out << __func__ << endl; double dExtents[6]; + double *origin, *spacing; + origin = sourceImage->Get(vtkDataObject::ORIGIN()); + spacing = sourceImage->Get(vtkDataObject::SPACING()); for(unsigned int i=0; i<6; i++) { // From source voxel coordinates to world coordinates - dExtents[i] = sourceImage->GetOrigin()[i/2] + sourceImage->GetSpacing()[i/2] * sourceExtent[i]; + dExtents[i] = origin[i/2] + spacing[i/2] * sourceExtent[i]; // From world coordinates to floating point target voxel coordinates dExtents[i] = (dExtents[i]- targetImage->GetOrigin()[i/2]) / targetImage->GetSpacing()[i/2]; diff --git a/vv/vvSlicer.h b/vv/vvSlicer.h index 1818e6c..43e0854 100644 --- a/vv/vvSlicer.h +++ b/vv/vvSlicer.h @@ -286,7 +286,7 @@ protected: private: void UpdateOrientation(); void UpdateDisplayExtent(); - void ConvertImageToImageDisplayExtent(vtkImageData *sourceImage, const int sourceExtent[6], + void ConvertImageToImageDisplayExtent(vtkInformation *sourceImage, const int sourceExtent[6], vtkImageData *targetImage, int targetExtent[6]); ///Sets the surfaces to be cut on the image slice: update the vtkCutter void SetContourSlice(); diff --git a/vv/vvToolRigidReg.cxx b/vv/vvToolRigidReg.cxx index 49b46f1..a3ee902 100644 --- a/vv/vvToolRigidReg.cxx +++ b/vv/vvToolRigidReg.cxx @@ -359,7 +359,6 @@ void vvToolRigidReg::SetTransform(vtkMatrix4x4 *matrix) vtkSmartPointer transform=vtkSmartPointer::New(); // TODO SR and BP: check on the list of transforms and not the first only mCurrentSlicerManager->GetImage()->GetTransform()[0]->SetMatrix(matrix); - //mCurrentSlicerManager->GetSlicer(2)->GetSlicingTransform()->SetMatrix(matrix); transform->Update(); Render(); dynamic_cast(mMainWindow)->ImageInfoChanged(); @@ -460,200 +459,8 @@ void vvToolRigidReg::ExtentMax(const double pointExtent[8][4], double maxExtent[ //------------------------------------------------------------------------------ void vvToolRigidReg::Render() -{ //out << __func__ << endl; -#if VTK_MAJOR_VERSION > 7 -double translationValues[4], translationValuesUpdate[4]; -mCurrentSlicerManager->GetImage()->GetTransform()[0]->Print(cout); -vtkMatrix4x4* matrix = mCurrentSlicerManager->GetImage()->GetTransform()[0]->GetMatrix(); -vtkMatrix4x4* matrixTranspose = matrix->NewInstance(); -for (int i=0; i<3; ++i) { - for (int j=0; j<3; ++j) - { - matrixTranspose->SetElement(i,j,matrix->GetElement(j,i)); - } -} -for (int j=0; j<3; ++j) { - translationValues[j] = matrix->GetElement(j,3); -} -translationValues[3] = 0.0; -matrix->MultiplyPoint(translationValues, translationValuesUpdate); -for (int i=0; i<3; ++i) { - matrixTranspose->SetElement(i,3,translationValuesUpdate[i]); -} -for (int i=0; i<4; ++i) { - matrixTranspose->SetElement(3,i,matrix->GetElement(3,i)); -} - -#endif for (int i=0; iGetNumberOfSlicers(); i++) { -#if VTK_MAJOR_VERSION > 7 - double pointExtent[8][4], pointExtentUpdate[8][4], pointOverlayExtent[8][4], pointOverlayExtentUpdate[8][4], centre[3], translation[3]; - std::vector w_ext; - w_ext=mCurrentSlicerManager->GetImage()->GetSize(); - pointExtent[0][0] = 0.0; - pointExtent[0][1] = 0.0; - pointExtent[0][2] = 0.0; - pointExtent[0][3] = 1.0; - pointExtent[1][0] = w_ext[0]-1; - pointExtent[1][1] = w_ext[1]-1; - pointExtent[1][2] = w_ext[2]-1; - pointExtent[1][3] = 1.0; - pointExtent[2][0] = 0.0; - pointExtent[2][1] = w_ext[1]-1; - pointExtent[2][2] = w_ext[2]-1; - pointExtent[2][3] = 1.0; - pointExtent[3][0] = w_ext[0]-1; - pointExtent[3][1] = 0.0; - pointExtent[3][2] = w_ext[2]-1; - pointExtent[3][3] = 1.0; - pointExtent[4][0] = w_ext[0]-1; - pointExtent[4][1] = w_ext[1]-1; - pointExtent[4][2] = 0.0; - pointExtent[4][3] = 1.0; - pointExtent[5][0] = 0.0; - pointExtent[5][1] = 0.0; - pointExtent[5][2] = w_ext[2]-1; - pointExtent[5][3] = 1.0; - pointExtent[6][0] = 0.0; - pointExtent[6][1] = w_ext[1]-1; - pointExtent[6][2] = 0.0; - pointExtent[6][3] = 1.0; - pointExtent[7][0] = w_ext[0]-1; - pointExtent[7][1] = 0.0; - pointExtent[7][2] = 0.0; - pointExtent[7][3] = 1.0; - - centre[0] = Xval->text().toDouble(); - centre[1] = Yval->text().toDouble(); - centre[2] = Zval->text().toDouble(); - - for (int k=0; k<8; ++k) { - for (int j=0; j<3; ++j) - { - pointOverlayExtent[k][j] = mCurrentSlicerManager->GetImage()->GetSpacing()[j]*pointExtent[k][j] - centre[j]; - } - pointOverlayExtent[k][3] = 0.0; - matrixTranspose->MultiplyPoint(pointOverlayExtent[k], pointOverlayExtentUpdate[k]); - for (int j=0; j<3; ++j) - { - pointOverlayExtentUpdate[k][j] = (pointOverlayExtentUpdate[k][j] + centre[j])/mCurrentSlicerManager->GetImage()->GetSpacing()[j]; - cout << pointOverlayExtentUpdate[k][j] << " "; - } - cout << endl; - } - cout << endl; - for (int k=0; k<8; ++k) { - for (int j=0; j<3; ++j) - { - pointExtent[k][j] = mCurrentSlicerManager->GetImage()->GetSpacing()[j] * pointExtent[k][j]; - } - matrixTranspose->MultiplyPoint(pointExtent[k], pointExtentUpdate[k]); - for (int j=0; j<3; ++j) - { - pointExtentUpdate[k][j] = (pointExtentUpdate[k][j])/mCurrentSlicerManager->GetImage()->GetSpacing()[j]; - cout << pointExtentUpdate[k][j] << " "; - } - cout << endl; - } - double extUpdateTemp[2][3], extOverlayUpdateTemp[2][3]; - int extUpdate[6]; - ExtentMax(pointExtentUpdate, extUpdateTemp); - ExtentMax(pointOverlayExtentUpdate, extOverlayUpdateTemp); - for (int j=0; j<3; ++j) { - extUpdate[2*j] = 0; - extUpdate[2*j+1] = itk::Math::Round(extUpdateTemp[1][j] - extUpdateTemp[0][j]); - } - mCurrentSlicerManager->GetSlicer(i)->SetRegisterExtent(extUpdate); - extUpdate[2*mCurrentSlicerManager->GetSlicer(i)->GetOrientation()] = mCurrentSlicerManager->GetSlicer(i)->GetSlice(); - extUpdate[2*mCurrentSlicerManager->GetSlicer(i)->GetOrientation()+1] = mCurrentSlicerManager->GetSlicer(i)->GetSlice(); - - vtkSmartPointer mapperOpenGL= vtkSmartPointer::New(); - try { - mapperOpenGL = dynamic_cast(mCurrentSlicerManager->GetSlicer(i)->GetImageActor()->GetMapper()); - } catch (const std::bad_cast& e) { - std::cerr << e.what() << std::endl; - std::cerr << "Conversion error" << std::endl; - return; - } - cout << extUpdate[0] << " " << extUpdate[1] << " " << extUpdate[2] << " " << extUpdate[3] << " " << extUpdate[4] << " " << extUpdate[5] << endl; - mapperOpenGL->SetCroppingRegion(extUpdate); - - if (mCurrentSlicerManager->GetSlicer(i)->GetOverlay() && mCurrentSlicerManager->GetSlicer(i)->GetOverlayActor()->GetVisibility()) { - int extOverlayUpdate[6]; - for (int j=0; j<3; ++j) { //Rotation - if (extOverlayUpdateTemp[1][j] - extOverlayUpdateTemp[0][j] > w_ext[j]-1) { - extOverlayUpdate[2*j] = 0; - extOverlayUpdate[2*j+1] = w_ext[j]-1; - } else { - extOverlayUpdate[2*j] = itk::Math::Round(extOverlayUpdateTemp[0][j]); - extOverlayUpdate[2*j+1] = itk::Math::Round(extOverlayUpdateTemp[1][j]); - } - } - - //Compute translation - double pointOrigin[4], pointOriginUpdate[4]; - for (int j=0; j<3; ++j) - { - pointOrigin[j] = 0 - centre[j]; - } - pointOrigin[3] = 0.0; - matrix->MultiplyPoint(pointOrigin, pointOriginUpdate); - for (int j=0; j<3; ++j) - { - pointOriginUpdate[j] = (pointOriginUpdate[j] + centre[j])/mCurrentSlicerManager->GetImage()->GetSpacing()[j]; - translation[j] = matrix->GetElement(j,3) - pointOriginUpdate[j]*mCurrentSlicerManager->GetImage()->GetSpacing()[j]; - pointOrigin[j] = translation[j]; - } - pointOrigin[3] = 0.0; - matrixTranspose->MultiplyPoint(pointOrigin, pointOriginUpdate); - for (int j=0; j<3; ++j) - { - translation[j] = pointOriginUpdate[j]/mCurrentSlicerManager->GetImage()->GetSpacing()[j]; - } - - for (int j=0; j<3; ++j) { //Translation - if (0 < extOverlayUpdateTemp[0][j] - translation[j]) { - extOverlayUpdate[2*j] = itk::Math::Round(extOverlayUpdateTemp[0][j] - translation[j]); - } else { - extOverlayUpdate[2*j] = 0; - } - if (extOverlayUpdateTemp[1][j] - translation[j] < w_ext[j]-1) { - extOverlayUpdate[2*j+1] = itk::Math::Round(extOverlayUpdateTemp[1][j] - translation[j]); - } else { - extOverlayUpdate[2*j+1] = w_ext[j]-1; - } - } - extOverlayUpdate[2*mCurrentSlicerManager->GetSlicer(i)->GetOrientation()] += mCurrentSlicerManager->GetSlicer(i)->GetSlice(); - extOverlayUpdate[2*mCurrentSlicerManager->GetSlicer(i)->GetOrientation()+1] = extOverlayUpdate[2*mCurrentSlicerManager->GetSlicer(i)->GetOrientation()]; - vtkSmartPointer mapperOpenGL= vtkSmartPointer::New(); - try { - mapperOpenGL = dynamic_cast(mCurrentSlicerManager->GetSlicer(i)->GetOverlayActor()->GetMapper()); - } catch (const std::bad_cast& e) { - std::cerr << e.what() << std::endl; - std::cerr << "Conversion error" << std::endl; - return; - } - - //Slice number - double spacing[4],spacingUpdate[4]; - int sliceNumber; - spacing[0] = 240*mCurrentSlicerManager->GetImage()->GetSpacing()[0]-centre[0]; - spacing[1] = 179*mCurrentSlicerManager->GetImage()->GetSpacing()[1]-centre[1]; - spacing[2] = 22*mCurrentSlicerManager->GetImage()->GetSpacing()[2]-centre[2]; - spacing[3] = 0; - matrixTranspose->MultiplyPoint(spacing, spacingUpdate); - spacingUpdate[0] = (spacingUpdate[0]+centre[0])/mCurrentSlicerManager->GetImage()->GetSpacing()[0]; - spacingUpdate[1] = (spacingUpdate[1]+centre[1])/mCurrentSlicerManager->GetImage()->GetSpacing()[1]; - spacingUpdate[2] = (spacingUpdate[2]+centre[2])/mCurrentSlicerManager->GetImage()->GetSpacing()[2]; - cout << spacingUpdate[0] << " " << spacingUpdate[1] << " " << spacingUpdate[2] << endl; - sliceNumber = mCurrentSlicerManager->GetSlicer(i)->GetSlice()*spacingUpdate[mCurrentSlicerManager->GetSlicer(i)->GetOrientation()]/mCurrentSlicerManager->GetImage()->GetSpacing()[mCurrentSlicerManager->GetSlicer(i)->GetOrientation()]; - extOverlayUpdate[2*mCurrentSlicerManager->GetSlicer(i)->GetOrientation()] = 12;//sliceNumber; - extOverlayUpdate[2*mCurrentSlicerManager->GetSlicer(i)->GetOrientation()+1] = 12;//sliceNumber; - - mapperOpenGL->SetCroppingRegion(extOverlayUpdate); - } -#endif mCurrentSlicerManager->GetSlicer(i)->ForceUpdateDisplayExtent(); mCurrentSlicerManager->GetSlicer(i)->Render(); }