fattori_g/reg23Topograms
0
0
Fork 0
Code Issues 1 Pull Requests Actions Packages Projects Releases Wiki Activity

R23 ROI Fixed and initial on LUT

Browse Source 
- Fixed bug in ROI itkRegion definition, missing division by spacing
- Consistent results using Small (1mm) or Large (2mm) FoV
- Some initial work on LUT, to polish...
This commit is contained in:
Proton local user
2023-05-19 23:40:46 +02:00
parent bc3d3e7e54
commit 47fe531aae
4 changed files with 76 additions and 102 deletions
Download Patch File Download Diff File Expand all files Collapse all files

14
reg23Topograms/itkDTRrecon/autoreg/itkNormalizedCorrelationTwoImageToOneImageMetric.hxx
View File

@ -188,15 +188,6 @@ NormalizedCorrelationTwoImageToOneImageMetric<TFixedImage, TMovingImage>::Calcul
++ti;
}
// sfm = sfm - m_meanCurr_f * sm - m_meanCurr_m * sf + this->m_NumberOfPixelsCounted * m_meanCurr_m * m_meanCurr_f;
// const RealType denom = -1 * sqrt(m_varianceCurr_f * m_varianceCurr_m) * this->m_NumberOfPixelsCounted;
// if (this->m_NumberOfPixelsCounted > 0 && denom != 0.0) {
// measure = sfm / denom;
// } else {
// measure = NumericTraits<MeasureType>::ZeroValue();
// }
if (this->m_SubtractMean && this->m_NumberOfPixelsCounted > 0) {
sff -= (sf * sf / this->m_NumberOfPixelsCounted);
smm -= (sm * sm / this->m_NumberOfPixelsCounted);
@ -273,6 +264,11 @@ NormalizedCorrelationTwoImageToOneImageMetric<TFixedImage, TMovingImage>::GetVal
// Calculate the measure value between fixed image 1 and the moving image
auto oldprecision = std::cout.precision();
// std::cout<<"Region " <<this->GetFixedImageRegion1() <<std::endl;
// std::cout<<"Buffered Region " <<this->GetFixedImage1()->GetBufferedRegion() <<std::endl;
// std::cout<<"Buffered Region " <<this->m_Filter1->GetOutput()->GetBufferedRegion() <<std::endl;
MeasureType measure1 = CalculateMeasure(1);
// std::cout.precision(std::numeric_limits<double>::digits10 + 2);
// std::cout << "Measure1: " << measure1 << " Counts: " << this->m_NumberOfPixelsCounted << std::endl;

6
reg23Topograms/itkDTRrecon/autoreg/itkgTwoImageToOneImageMetric.hxx
View File

@ -250,6 +250,12 @@ void gTwoImageToOneImageMetric<TFixedImage, TMovingImage>::Initialize()
m_FixedImage2->GetSource()->Update();
}
//std::cout<<"FixedImageRegion1 " <<this->GetFixedImageRegion1() <<std::endl;
//std::cout<<"m_FixedImage1 buffered " <<this->m_FixedImage1->GetBufferedRegion() <<std::endl;
//std::cout<<"Moving buffered " <<this->m_Filter1->GetOutput()->GetBufferedRegion() <<std::endl;
//std::cout<<"GetFixedImage1 buffered " <<this->GetFixedImage1()->GetBufferedRegion() <<std::endl;
//std::cout<<"Filter1 buffered " <<this->m_Filter1->GetOutput()->GetBufferedRegion() <<std::endl;
// Make sure the FixedImageRegion is within the FixedImage buffered region
if (!m_FixedImageRegion1.Crop(m_FixedImage1->GetBufferedRegion())) {
itkExceptionMacro(<< "FixedImageRegion1 does not overlap the fixed image buffered region");

150
reg23Topograms/itkDTRrecon/itkImageProcessor.cpp
View File

@ -1968,44 +1968,45 @@ vtkImageData* itkImageProcessor::GetLocalizer1VTK()
toVTKLocalizer1->SetInput(rescaler1->GetOutput());
toVTKLocalizer1->Update();
if(false) {
using ImageRegionType3D = ImageType3D::RegionType;
using SizeType3D = ImageRegionType3D::SizeType;
using ImageDirectionType3D = ImageType3D::DirectionType;
// if(false) {
//// using ImageRegionType3D = ImageType3D::RegionType;
//// using SizeType3D = ImageRegionType3D::SizeType;
//// using ImageDirectionType3D = ImageType3D::DirectionType;
ImageRegionType3D region3D = rescaler1->GetOutput()->GetBufferedRegion();
SizeType3D size3D = region3D.GetSize();
ImageType3D::PointType origin3D = rescaler1->GetOutput()->GetOrigin();
const itk::Vector<double, 3> resolution3D = rescaler1->GetOutput()->GetSpacing();
ImageDirectionType3D imagDirection = rescaler1->GetOutput()->GetDirection();
//// ImageRegionType3D region3D = rescaler1->GetOutput()->GetBufferedRegion();
//// SizeType3D size3D = region3D.GetSize();
//// //ImageType3D::PointType origin3D = rescaler1->GetOutput()->GetOrigin();
//// const itk::Vector<double, 3> resolution3D = rescaler1->GetOutput()->GetSpacing();
//// ImageDirectionType3D imagDirection = rescaler1->GetOutput()->GetDirection();
/* calculate image size in 3D Space */
using VectorType = itk::Vector<double, 3>;
VectorType Dest;
Dest[0]=(size3D[0]-1) * resolution3D[0];
Dest[1]=(size3D[1]-1) * resolution3D[1];
Dest[2]=(size3D[2]-1) * resolution3D[2];
//// /* calculate image size in 3D Space */
//// using VectorType = itk::Vector<double, 3>;
//// VectorType Dest;
//// Dest[0]=(size3D[0]-1) * resolution3D[0];
//// Dest[1]=(size3D[1]-1) * resolution3D[1];
//// Dest[2]=(size3D[2]-1) * resolution3D[2];
ImageType3D::PointType LastVoxelPosition =
origin3D + (imagDirection * Dest);
// std::cout<<" -------- Localizer 1 ITK --------"<<std::endl;
// std::cout<<"size3D " <<size3D <<std::endl;
// std::cout<<"resolution3D " <<resolution3D <<std::endl;
// std::cout<<"imagDirection " <<imagDirection <<std::endl;
// std::cout<<"First Voxel position " <<origin3D <<std::endl;
// std::cout<<"Last Voxel position: "<<LastVoxelPosition<<std::endl;
// // ImageType3D::PointType LastVoxelPosition =
// // origin3D + (imagDirection * Dest);
// // std::cout<<" -------- Localizer 1 ITK --------"<<std::endl;
// // std::cout<<"size3D " <<size3D <<std::endl;
// // std::cout<<"resolution3D " <<resolution3D <<std::endl;
// // std::cout<<"imagDirection " <<imagDirection <<std::endl;
// // std::cout<<"First Voxel position " <<origin3D <<std::endl;
// // std::cout<<"Last Voxel position: "<<LastVoxelPosition<<std::endl;
double* dBounds = toVTKLocalizer1->GetOutput()->GetBounds();
// // double* dBounds = toVTKLocalizer1->GetOutput()->GetBounds();
// std::cout<< "-------- Localizer 1 VTK --------" <<std::endl;
// std::cout<<"Bounds: "<<dBounds[0]<<" "<<dBounds[1]<<" "
// <<dBounds[2]<<" "<<dBounds[3]<<" "
// <<dBounds[4]<<" "<<dBounds[5]<<std::endl;
// std::cout<< "-------- -------- --------" <<std::endl;
}
// // std::cout<< "-------- Localizer 1 VTK --------" <<std::endl;
// // std::cout<<"Bounds: "<<dBounds[0]<<" "<<dBounds[1]<<" "
// // <<dBounds[2]<<" "<<dBounds[3]<<" "
// // <<dBounds[4]<<" "<<dBounds[5]<<std::endl;
// // std::cout<< "-------- -------- --------" <<std::endl;
// }
return
toVTKLocalizer1->GetOutput();
@ -2030,42 +2031,42 @@ vtkImageData* itkImageProcessor::GetLocalizer2VTK()
if(true) {
using ImageRegionType3D = ImageType3D::RegionType;
using SizeType3D = ImageRegionType3D::SizeType;
using ImageDirectionType3D = ImageType3D::DirectionType;
// if(true) {
// using ImageRegionType3D = ImageType3D::RegionType;
// using SizeType3D = ImageRegionType3D::SizeType;
// using ImageDirectionType3D = ImageType3D::DirectionType;
ImageRegionType3D region3D = rescaler2->GetOutput()->GetBufferedRegion();
SizeType3D size3D = region3D.GetSize();
ImageType3D::PointType origin3D = rescaler2->GetOutput()->GetOrigin();
const itk::Vector<double, 3> resolution3D = rescaler2->GetOutput()->GetSpacing();
ImageDirectionType3D imagDirection = rescaler2->GetOutput()->GetDirection();
// ImageRegionType3D region3D = rescaler2->GetOutput()->GetBufferedRegion();
// SizeType3D size3D = region3D.GetSize();
// ImageType3D::PointType origin3D = rescaler2->GetOutput()->GetOrigin();
// const itk::Vector<double, 3> resolution3D = rescaler2->GetOutput()->GetSpacing();
// ImageDirectionType3D imagDirection = rescaler2->GetOutput()->GetDirection();
/* calculate image size in 3D Space */
using VectorType = itk::Vector<double, 3>;
VectorType Dest;
Dest[0]=(size3D[0]-1) * resolution3D[0];
Dest[1]=(size3D[1]-1) * resolution3D[1];
Dest[2]=(size3D[2]-1) * resolution3D[2];
// /* calculate image size in 3D Space */
// using VectorType = itk::Vector<double, 3>;
// VectorType Dest;
// Dest[0]=(size3D[0]-1) * resolution3D[0];
// Dest[1]=(size3D[1]-1) * resolution3D[1];
// Dest[2]=(size3D[2]-1) * resolution3D[2];
ImageType3D::PointType LastVoxelPosition =
origin3D + (imagDirection * Dest);
// ImageType3D::PointType LastVoxelPosition =
// origin3D + (imagDirection * Dest);
// std::cout<<" -------- Localizer 2 ITK --------"<<std::endl;
// std::cout<<"size3D " <<size3D <<std::endl;
// std::cout<<"resolution3D " <<resolution3D <<std::endl;
// std::cout<<"imagDirection " <<imagDirection <<std::endl;
// std::cout<<"First Voxel position " <<origin3D <<std::endl;
// std::cout<<"Last Voxel position: "<<LastVoxelPosition<<std::endl;
// // std::cout<<" -------- Localizer 2 ITK --------"<<std::endl;
// // std::cout<<"size3D " <<size3D <<std::endl;
// // std::cout<<"resolution3D " <<resolution3D <<std::endl;
// // std::cout<<"imagDirection " <<imagDirection <<std::endl;
// // std::cout<<"First Voxel position " <<origin3D <<std::endl;
// // std::cout<<"Last Voxel position: "<<LastVoxelPosition<<std::endl;
// double* dBounds = toVTKLocalizer2->GetOutput()->GetBounds();
// std::cout<< "-------- Localizer 2 VTK --------" <<std::endl;
// std::cout<<"Bounds: "<<dBounds[0]<<" "<<dBounds[1]<<" "
// <<dBounds[2]<<" "<<dBounds[3]<<" "
// <<dBounds[4]<<" "<<dBounds[5]<<std::endl;
// std::cout<< "-------- -------- --------" <<std::endl;
// // double* dBounds = toVTKLocalizer2->GetOutput()->GetBounds();
// // std::cout<< "-------- Localizer 2 VTK --------" <<std::endl;
// // std::cout<<"Bounds: "<<dBounds[0]<<" "<<dBounds[1]<<" "
// // <<dBounds[2]<<" "<<dBounds[3]<<" "
// // <<dBounds[4]<<" "<<dBounds[5]<<std::endl;
// // std::cout<< "-------- -------- --------" <<std::endl;
}
// }
@ -2154,11 +2155,6 @@ vtkImageData* itkImageProcessor::GetProjection1VTK()
vtkMatrix4x4 * itkImageProcessor::GetProjection1VTKTransform()
{
// std::cout<< "Composed " <<std::endl;
// interpolator1->GetComposedTransform()->Print(std::cout);
// std::cout<< "Inverse " <<std::endl;
// interpolator1->GetInverseTransform()->Print(std::cout);
m_Projection1VTKTransform->Identity();
for(int iIdx = 0; iIdx<3 ; iIdx++){
@ -2398,22 +2394,6 @@ itkImageProcessor::GetCompleteIsocentricTransform
return nullptr;
}
// std::cout<<" --- GetCompleteIsocentricTransform --- "<<std::endl;
// std::cout<< "R: "<<
// m_TransformMetaInfo->GetR()[0] <<" "<<
// m_TransformMetaInfo->GetR()[1] <<" "<<
// m_TransformMetaInfo->GetR()[2] <<" "<<std::endl;
// std::cout<< "T: "<<
// m_TransformMetaInfo->GetT()[0] <<" "<<
// m_TransformMetaInfo->GetT()[1] <<" "<<
// m_TransformMetaInfo->GetT()[2] <<" "<<std::endl;
// std::cout<< "Import Offset: "<<
// m_CTMetaInfo->GetImportOffset()[0]<<" "<<
// m_CTMetaInfo->GetImportOffset()[1]<<" "<<
// m_CTMetaInfo->GetImportOffset()[2]<<" "<<std::endl;
// std::cout<<" --- --- --- "<<std::endl;
ImageType3D::PointType mISORTIEC =
m_CTMetaInfo->GetLPS2IECDirections() * m_RTMetaInfo->GetIsocenterLPS();
@ -2465,8 +2445,8 @@ void itkImageProcessor::SetRegionFixed1(
roiAutoReg1.SetIndex(index1);
roiAutoReg1.SetSize(size1);
std::cout << "itkImageProcessor " << std::endl;
std::cout << roiAutoReg1 << std::endl;
//std::cout << "itkImageProcessor " << std::endl;
//std::cout << roiAutoReg1 << std::endl;
this->m_R23->SetroiAutoReg1(roiAutoReg1);
@ -2490,8 +2470,8 @@ void itkImageProcessor::SetRegionFixed2(
roiAutoReg2.SetIndex(index2);
roiAutoReg2.SetSize(size2);
std::cout << "itkImageProcessor " << std::endl;
std::cout << roiAutoReg2 << std::endl;
//std::cout << "itkImageProcessor " << std::endl;
//std::cout << roiAutoReg2 << std::endl;
this->m_R23->SetroiAutoReg2(roiAutoReg2);

8
reg23Topograms/itkDTRrecon/itkReg23.cpp
View File

@ -210,9 +210,6 @@ void itkReg23::InitializeRegistration(
registration->SetFixedImageRegion1(m_roiAutoReg1);
registration->SetFixedImageRegion2(m_roiAutoReg2);
std::cout<< "Using Region1:"<<m_roiAutoReg1 <<std::endl;
std::cout<< "Using Region2:"<<m_roiAutoReg2 <<std::endl;
registration->SetTransformMetaInfo(m_r23Meta);
@ -391,11 +388,6 @@ double itkReg23::GetCurrentMetricValue()
}
//void itkReg23::SetFullROI(bool fullROI)
//{
// m_UseFullROI = fullROI;
// // TODO: Remove this function when ROI functinalitz has been implemented.
//}
}