Restored TestApp functionality with simplifed rendering pipeline
This commit is contained in:
@ -84,24 +84,11 @@ const char *gGetStringValueFromTag(const gdcm::Tag& t, const gdcm::DataSet& ds)
|
||||
itkImageProcessor::itkImageProcessor()
|
||||
{
|
||||
|
||||
std::cout<<"itkImageProcessor::itkImageProcessor contructor"<<std::endl;
|
||||
metric = MetricType::New();
|
||||
// std::cout<<"itkImageProcessor::itkImageProcessor contructor"<<std::endl;
|
||||
transform = TransformType::New();
|
||||
optimizer = OptimizerType::New();
|
||||
interpolator1 = InterpolatorType::New();
|
||||
interpolator2 = InterpolatorType::New();
|
||||
registration = RegistrationType::New();
|
||||
|
||||
metric->ComputeGradientOff();
|
||||
metric->SetSubtractMean(true);
|
||||
|
||||
registration->SetMetric(metric);
|
||||
registration->SetOptimizer(optimizer);
|
||||
registration->SetTransform(transform);
|
||||
registration->SetInterpolator1(interpolator1);
|
||||
registration->SetInterpolator2(interpolator2);
|
||||
|
||||
//imageReader2D = ImageReaderType2D::New();
|
||||
imageReader3D = ImageReaderType3D::New();
|
||||
imageSeriesReader3D = ImageSeriesReaderType3D::New();
|
||||
|
||||
@ -132,25 +119,6 @@ itkImageProcessor::itkImageProcessor()
|
||||
image2D2Loaded= false;
|
||||
image3DLoaded= false;
|
||||
|
||||
flipFilter1 = FlipFilterType::New();
|
||||
flipFilter2 = FlipFilterType::New();
|
||||
flipFilterLocalizer1 = FlipFilterType::New();
|
||||
flipFilterLocalizer2 = FlipFilterType::New();
|
||||
|
||||
using FlipAxesArrayType = FlipFilterType::FlipAxesArrayType;
|
||||
FlipAxesArrayType flipArray;
|
||||
flipArray[0] = false;
|
||||
flipArray[1] = false;
|
||||
flipArray[2] = false;
|
||||
// flipArray[0] = true;
|
||||
// flipArray[1] = true;
|
||||
// flipArray[2] = false;
|
||||
flipFilter1->SetFlipAxes(flipArray);
|
||||
// flipArray[0] = false;
|
||||
// flipArray[1] = true;
|
||||
// flipArray[2] = false;
|
||||
flipFilter2->SetFlipAxes(flipArray);
|
||||
|
||||
caster3D = CastFilterType3D::New();
|
||||
caster2D1 = CastFilterType3D::New();
|
||||
caster2D2 = CastFilterType3D::New();
|
||||
@ -160,8 +128,6 @@ itkImageProcessor::itkImageProcessor()
|
||||
toVTKLocalizer1 = ITKtoVTKFilterType::New();
|
||||
toVTKLocalizer2 = ITKtoVTKFilterType::New();
|
||||
|
||||
|
||||
|
||||
m_3DPatOrientation = eImageOrientationType::NotDefined;
|
||||
|
||||
m_LATSourceDupli = DuplicatorType::New();
|
||||
@ -236,7 +202,7 @@ void itkImageProcessor::SetCustom_2Dcenter(double dX1,double dY1, double dX2,dou
|
||||
int itkImageProcessor::load3DSerie(const char * pcDirName)
|
||||
{
|
||||
|
||||
std::cout<<"itkImageProcessor::load3DSerie"<<std::endl;
|
||||
// std::cout<<"itkImageProcessor::load3DSerie"<<std::endl;
|
||||
|
||||
using NamesGeneratorType = itk::GDCMSeriesFileNames;
|
||||
NamesGeneratorType::Pointer nameGenerator = NamesGeneratorType::New();
|
||||
@ -443,7 +409,7 @@ double
|
||||
tPatOrientation pOrient,
|
||||
double pAngleIEC){
|
||||
|
||||
// std::cout<< "CalcProjectionAngleLPS :: pAngleIEC " << pAngleIEC <<std::endl;
|
||||
// std::cout<< "CalcProjectèionAngleLPS :: pAngleIEC " << pAngleIEC <<std::endl;
|
||||
|
||||
double dProj = pAngleIEC;
|
||||
|
||||
@ -710,13 +676,7 @@ int itkImageProcessor::load2D(const char * pcFName){
|
||||
|
||||
const InternalImageType::DirectionType newDirection =
|
||||
MyLocalizerImage->GetDirection() * image2D1InDirectionInv;
|
||||
// const InternalImageType::DirectionType newDirection =
|
||||
// toIECDirection * MyLocalizerImage->GetDirection() ;
|
||||
|
||||
//filter->SetOutputDirection(newDirection);
|
||||
//filter->ChangeDirectionOn();
|
||||
//filter->SetOutputOrigin(image2Dcenter_3DCOV);
|
||||
//filter->ChangeOriginOn();
|
||||
try
|
||||
{
|
||||
filter->UpdateOutputInformation();
|
||||
@ -730,7 +690,6 @@ int itkImageProcessor::load2D(const char * pcFName){
|
||||
filter->Update();
|
||||
|
||||
|
||||
InternalImageType::Pointer m_Image;
|
||||
double* m_ImageIntensity;
|
||||
FlipFilterType::Pointer m_FlipFilter;
|
||||
DuplicatorType::Pointer m_Duplicator;
|
||||
@ -739,20 +698,16 @@ int itkImageProcessor::load2D(const char * pcFName){
|
||||
switch (currProjOrientation) {
|
||||
case eProjectionOrientationType::PA:
|
||||
|
||||
m_Image = image2D1In.GetPointer();
|
||||
m_ImageIntensity = image1IntensityWindow;
|
||||
m_FlipFilter = flipFilterLocalizer1.GetPointer();
|
||||
m_Duplicator = m_PASourceDupli.GetPointer();
|
||||
m_Duplicator = m_PASourceDupli;
|
||||
m_ImLoaded = &image2D1Loaded;
|
||||
|
||||
break;
|
||||
|
||||
case eProjectionOrientationType::LAT:
|
||||
|
||||
m_Image = image2D2In.GetPointer();
|
||||
m_ImageIntensity = image2IntensityWindow;
|
||||
m_FlipFilter = flipFilterLocalizer2.GetPointer();
|
||||
m_Duplicator = m_LATSourceDupli.GetPointer();
|
||||
m_Duplicator = m_LATSourceDupli;
|
||||
m_ImLoaded = &image2D2Loaded;
|
||||
|
||||
break;
|
||||
@ -767,9 +722,8 @@ int itkImageProcessor::load2D(const char * pcFName){
|
||||
|
||||
m_Duplicator->SetInputImage(filter->GetOutput());
|
||||
m_Duplicator->Update();
|
||||
m_Image = m_Duplicator->GetOutput();
|
||||
|
||||
m_FlipFilter->SetInput(m_Image);
|
||||
// m_FlipFilter->SetInput(m_Image);
|
||||
|
||||
|
||||
return
|
||||
@ -833,7 +787,7 @@ void itkImageProcessor::InitializeProjector()
|
||||
isocenter3D[0] = origin3D[0] + resolution3D[0] * isocenter[0];
|
||||
isocenter3D[1] = origin3D[1] + resolution3D[1] * isocenter[1];
|
||||
isocenter3D[2] = origin3D[2] + resolution3D[2] * isocenter[2];
|
||||
std::cout<<"Isocenter location given by the user"<<std::endl;
|
||||
// std::cout<<"Isocenter location given by the user"<<std::endl;
|
||||
}
|
||||
else
|
||||
{
|
||||
@ -843,10 +797,10 @@ void itkImageProcessor::InitializeProjector()
|
||||
isocenter3D[2] = origin3D[2] + resolution3D[2] * static_cast<double>(size3D[2]) / 2.0;
|
||||
}
|
||||
|
||||
std::cout<<" -------- PROJECTOR GEOMETRY --------"<<std::endl;
|
||||
std::cout<<"Projector Isocenter "<<isocenter3D[0]<<" "<<isocenter3D[1]<<" "<<isocenter3D[2] <<std::endl;
|
||||
std::cout<<"Projector scd "<<scd <<std::endl;
|
||||
std::cout<<" -------- -------- --------"<<std::endl;
|
||||
// std::cout<<" -------- PROJECTOR GEOMETRY --------"<<std::endl;
|
||||
// std::cout<<"Projector Isocenter "<<isocenter3D[0]<<" "<<isocenter3D[1]<<" "<<isocenter3D[2] <<std::endl;
|
||||
// std::cout<<"Projector scd "<<scd <<std::endl;
|
||||
// std::cout<<" -------- -------- --------"<<std::endl;
|
||||
|
||||
transform->SetCenter(isocenter3D);
|
||||
// 2D Image 1
|
||||
@ -981,9 +935,9 @@ void itkImageProcessor::GetProjectionImages(){
|
||||
|
||||
if(false){
|
||||
// if(image2D1Loaded){
|
||||
resampleFilter1->SetSize(image2D1In->GetLargestPossibleRegion().GetSize());
|
||||
resampleFilter1->SetOutputOrigin(image2D1In->GetOrigin());
|
||||
resampleFilter1->SetOutputSpacing(image2D1In->GetSpacing());
|
||||
// resampleFilter1->SetSize(m_Input2DImage_PA ->GetLargestPossibleRegion().GetSize());
|
||||
// resampleFilter1->SetOutputOrigin(m_Input2DImage_PA ->GetOrigin());
|
||||
// resampleFilter1->SetOutputSpacing(m_Input2DImage_PA ->GetSpacing());
|
||||
|
||||
// std::cout<<"Size "<<caster2D1->GetOutput()->GetLargestPossibleRegion().GetSize()<<std::endl;
|
||||
// std::cout<<"Origin "<<caster2D1->GetOutput()->GetOrigin()<<std::endl;
|
||||
@ -1025,9 +979,9 @@ void itkImageProcessor::GetProjectionImages(){
|
||||
|
||||
if(false){
|
||||
// if(image2D2Loaded){
|
||||
resampleFilter2->SetSize(image2D2In->GetLargestPossibleRegion().GetSize());
|
||||
resampleFilter2->SetOutputOrigin(image2D2In->GetOrigin());
|
||||
resampleFilter2->SetOutputSpacing(image2D2In->GetSpacing());
|
||||
// resampleFilter2->SetSize(m_Input2DImage_LAT->GetLargestPossibleRegion().GetSize());
|
||||
// resampleFilter2->SetOutputOrigin(m_Input2DImage_LAT ->GetOrigin());
|
||||
// resampleFilter2->SetOutputSpacing(m_Input2DImage_LAT ->GetSpacing());
|
||||
|
||||
} else {
|
||||
|
||||
@ -1122,13 +1076,8 @@ void itkImageProcessor::GetProjectionImages(){
|
||||
filter1->Update();
|
||||
filter2->Update();
|
||||
|
||||
// As explained before, the computed projection is upsided-down.
|
||||
// Here we use a FilpImageFilter to flip the images in y-direction.
|
||||
/* no flip required because now we handle properly the image directions. */
|
||||
flipFilter1->SetInput(filter1->GetOutput());
|
||||
flipFilter2->SetInput(filter2->GetOutput());
|
||||
flipFilter1->Update();
|
||||
flipFilter2->Update();
|
||||
imageDRT1In = filter1->GetOutput();
|
||||
imageDRT2In = filter2->GetOutput();
|
||||
|
||||
}
|
||||
|
||||
@ -1146,7 +1095,7 @@ void itkImageProcessor::Write2DImages(){
|
||||
RescaleFilterType::Pointer rescaler1 = RescaleFilterType::New();
|
||||
rescaler1->SetOutputMinimum(0);
|
||||
rescaler1->SetOutputMaximum(255);
|
||||
rescaler1->SetInput(image2D1In);
|
||||
rescaler1->SetInput(m_PASourceDupli->GetOutput());
|
||||
|
||||
WriterType::Pointer writer1 = WriterType::New();
|
||||
|
||||
@ -1171,7 +1120,7 @@ void itkImageProcessor::Write2DImages(){
|
||||
RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New();
|
||||
rescaler2->SetOutputMinimum(0);
|
||||
rescaler2->SetOutputMaximum(255);
|
||||
rescaler2->SetInput(image2D2In);
|
||||
rescaler2->SetInput(m_LATSourceDupli->GetOutput());
|
||||
|
||||
WriterType::Pointer writer2 = WriterType::New();
|
||||
|
||||
@ -1200,7 +1149,7 @@ vtkImageData* itkImageProcessor::GetLocalizer1VTK()
|
||||
RescaleFilterType::Pointer rescaler1 = RescaleFilterType::New();
|
||||
rescaler1->SetOutputMinimum(0);
|
||||
rescaler1->SetOutputMaximum(255);
|
||||
rescaler1->SetInput(flipFilterLocalizer1->GetOutput());
|
||||
rescaler1->SetInput(m_PASourceDupli->GetOutput());
|
||||
|
||||
rescaler1->Update();
|
||||
|
||||
@ -1252,13 +1201,15 @@ vtkImageData* itkImageProcessor::GetLocalizer1VTK()
|
||||
|
||||
vtkImageData* itkImageProcessor::GetLocalizer2VTK()
|
||||
{
|
||||
|
||||
|
||||
// Rescale the intensity of the projection images to 0-255 for output.
|
||||
using RescaleFilterType = itk::RescaleIntensityImageFilter<InternalImageType, OutputImageType>;
|
||||
|
||||
RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New();
|
||||
rescaler2->SetOutputMinimum(0);
|
||||
rescaler2->SetOutputMaximum(255);
|
||||
rescaler2->SetInput(flipFilterLocalizer2->GetOutput());
|
||||
rescaler2->SetInput(m_LATSourceDupli ->GetOutput());
|
||||
|
||||
rescaler2->Update();
|
||||
|
||||
@ -1288,19 +1239,19 @@ vtkImageData* itkImageProcessor::GetLocalizer2VTK()
|
||||
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;
|
||||
|
||||
}
|
||||
|
||||
@ -1318,7 +1269,8 @@ vtkImageData* itkImageProcessor::GetProjection1VTK()
|
||||
RescaleFilterType::Pointer rescaler1 = RescaleFilterType::New();
|
||||
rescaler1->SetOutputMinimum(0);
|
||||
rescaler1->SetOutputMaximum(255);
|
||||
rescaler1->SetInput(flipFilter1->GetOutput());
|
||||
rescaler1->SetInput( imageDRT1In );
|
||||
// flipFilter1->GetOutput());
|
||||
|
||||
rescaler1->Update();
|
||||
|
||||
@ -1326,41 +1278,41 @@ vtkImageData* itkImageProcessor::GetProjection1VTK()
|
||||
toVTK2D1->Update();
|
||||
|
||||
|
||||
using ImageRegionType3D = ImageType3D::RegionType;
|
||||
using SizeType3D = ImageRegionType3D::SizeType;
|
||||
using ImageDirectionType3D = ImageType3D::DirectionType;
|
||||
// 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);
|
||||
// ImageType3D::PointType LastVoxelPosition =
|
||||
// origin3D + (imagDirection * Dest);
|
||||
|
||||
std::cout<<" -------- Projection 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;
|
||||
// std::cout<<" -------- Projection 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 = toVTK2D1->GetOutput()->GetBounds();
|
||||
// double* dBounds = toVTK2D1->GetOutput()->GetBounds();
|
||||
|
||||
std::cout<< "-------- Proj 1 --------" <<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<< "-------- Proj 1 --------" <<std::endl;
|
||||
// std::cout<<"Bounds: "<<dBounds[0]<<" "<<dBounds[1]<<" "
|
||||
// <<dBounds[2]<<" "<<dBounds[3]<<" "
|
||||
// <<dBounds[4]<<" "<<dBounds[5]<<std::endl;
|
||||
// std::cout<< "-------- -------- --------" <<std::endl;
|
||||
|
||||
return
|
||||
toVTK2D1->GetOutput();
|
||||
@ -1374,50 +1326,51 @@ vtkImageData* itkImageProcessor::GetProjection2VTK()
|
||||
RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New();
|
||||
rescaler2->SetOutputMinimum(0);
|
||||
rescaler2->SetOutputMaximum(255);
|
||||
rescaler2->SetInput(flipFilter2->GetOutput());
|
||||
rescaler2->SetInput( imageDRT2In );
|
||||
// flipFilter2->GetOutput());
|
||||
|
||||
rescaler2->Update();
|
||||
|
||||
rescaler2->Print(std::cout);
|
||||
// rescaler2->Print(std::cout);
|
||||
|
||||
using ImageRegionType3D = ImageType3D::RegionType;
|
||||
using SizeType3D = ImageRegionType3D::SizeType;
|
||||
using ImageDirectionType3D = ImageType3D::DirectionType;
|
||||
// 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];
|
||||
// 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<<" -------- Projection 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<<" -------- Projection 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;
|
||||
|
||||
|
||||
toVTK2D2->SetInput(rescaler2->GetOutput());
|
||||
toVTK2D2->Update();
|
||||
|
||||
double* dBounds = toVTK2D2->GetOutput()->GetBounds();
|
||||
// double* dBounds = toVTK2D2->GetOutput()->GetBounds();
|
||||
|
||||
std::cout<< "-------- Proj 2 --------" <<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<< "-------- Proj 2 --------" <<std::endl;
|
||||
// std::cout<<"Bounds: "<<dBounds[0]<<" "<<dBounds[1]<<" "
|
||||
// <<dBounds[2]<<" "<<dBounds[3]<<" "
|
||||
// <<dBounds[4]<<" "<<dBounds[5]<<std::endl;
|
||||
// std::cout<< "-------- -------- --------" <<std::endl;
|
||||
|
||||
return
|
||||
toVTK2D2->GetOutput();
|
||||
@ -1432,12 +1385,14 @@ void itkImageProcessor::WriteProjectionImages()
|
||||
RescaleFilterType::Pointer rescaler1 = RescaleFilterType::New();
|
||||
rescaler1->SetOutputMinimum(0);
|
||||
rescaler1->SetOutputMaximum(255);
|
||||
rescaler1->SetInput(flipFilter1->GetOutput());
|
||||
rescaler1->SetInput( imageDRT1In );
|
||||
// flipFilter1->GetOutput());
|
||||
|
||||
RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New();
|
||||
rescaler2->SetOutputMinimum(0);
|
||||
rescaler2->SetOutputMaximum(255);
|
||||
rescaler2->SetInput(flipFilter2->GetOutput());
|
||||
rescaler2->SetInput( imageDRT2In );
|
||||
// flipFilter2->GetOutput());
|
||||
|
||||
|
||||
using WriterType = itk::ImageFileWriter<OutputImageType>;
|
||||
@ -1475,486 +1430,6 @@ void itkImageProcessor::WriteProjectionImages()
|
||||
}
|
||||
|
||||
|
||||
//int itkImageProcessor::InitializeRegistration(){
|
||||
|
||||
|
||||
// // image2D1Rescaler->Update();
|
||||
// // image2D2Rescaler->Update();
|
||||
|
||||
// using FlipAxesArrayType = FlipFilterType::FlipAxesArrayType;
|
||||
|
||||
// FlipAxesArrayType flipArray;
|
||||
// flipArray[0] = true;
|
||||
// flipArray[1] = true;
|
||||
// flipArray[2] = false;
|
||||
// flipFilterLocalizer1->SetFlipAxes(flipArray);
|
||||
// flipFilterLocalizer1->Update();
|
||||
// flipArray[0] = false;
|
||||
// flipArray[1] = true;
|
||||
// flipArray[2] = false;
|
||||
// flipFilterLocalizer2->SetFlipAxes(flipArray);
|
||||
// flipFilterLocalizer2->Update();
|
||||
|
||||
// caster3D->Update();
|
||||
// using RescaleFilterType = itk::RescaleIntensityImageFilter<InternalImageType, InternalImageType>;
|
||||
|
||||
// RescaleFilterType::Pointer rescaler1 = RescaleFilterType::New();
|
||||
// rescaler1->SetOutputMinimum(0);
|
||||
// rescaler1->SetOutputMaximum(255);
|
||||
// rescaler1->SetInput(flipFilterLocalizer1->GetOutput());
|
||||
|
||||
// rescaler1->Update();
|
||||
|
||||
// using RescaleFilterType = itk::RescaleIntensityImageFilter<InternalImageType, InternalImageType>;
|
||||
|
||||
// RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New();
|
||||
// rescaler2->SetOutputMinimum(0);
|
||||
// rescaler2->SetOutputMaximum(255);
|
||||
// rescaler2->SetInput(flipFilterLocalizer2->GetOutput());
|
||||
|
||||
// rescaler2->Update();
|
||||
|
||||
|
||||
|
||||
// // metric->SetFixedImage1(rescaler1->GetOutput());
|
||||
// // metric->SetFixedImage2(rescaler2->GetOutput());
|
||||
// // metric->SetMovingImage(caster3D->GetOutput());
|
||||
// // metric->SetInterpolator1(interpolator1);
|
||||
// // metric->SetInterpolator2(interpolator2);
|
||||
|
||||
|
||||
// // TransformType::OutputVectorType translation;
|
||||
|
||||
// // translation[0] = 0.;
|
||||
// // translation[1] = 0.;
|
||||
// // translation[2] = 0.;
|
||||
|
||||
// // transform->SetTranslation(translation);
|
||||
|
||||
// // // constant for converting degrees to radians
|
||||
// // const double dtr = (atan(1.0) * 4.0) / 180.0;
|
||||
// // transform->SetRotation(dtr * 0., dtr * 0., dtr * 0.);
|
||||
|
||||
// // metric->SetTransform(transform);
|
||||
|
||||
// // std::cout<<" Curr VALUE: " <<
|
||||
// // metric->GetValue(transform->GetParameters())
|
||||
// // <<std::endl;
|
||||
// // return 0;
|
||||
|
||||
|
||||
// std::cout<<
|
||||
// "2d1 "<<
|
||||
// rescaler1->GetOutput()->GetDirection()
|
||||
// <<std::endl;
|
||||
// std::cout<<
|
||||
// "2d2 "<<
|
||||
// rescaler2->GetOutput()->GetDirection()
|
||||
// <<std::endl;
|
||||
// std::cout<<
|
||||
// "3d "<<
|
||||
// caster3D->GetOutput()->GetDirection()
|
||||
// <<std::endl;
|
||||
|
||||
|
||||
// registration->SetFixedImage1(rescaler1->GetOutput());
|
||||
// registration->SetFixedImage2(rescaler2->GetOutput());
|
||||
// registration->SetMovingImage(caster3D->GetOutput());
|
||||
|
||||
// // Initialise the transform
|
||||
// // ~~~~~~~~~~~~~~~~~~~~~~~~
|
||||
|
||||
// // Set the order of the computation. Default ZXY
|
||||
// transform->SetComputeZYX(true);
|
||||
|
||||
|
||||
// // The transform is initialised with the translation [tx,ty,tz] and
|
||||
// // rotation [rx,ry,rz] specified on the command line
|
||||
|
||||
// TransformType::OutputVectorType translation;
|
||||
|
||||
// translation[0] = 0.;
|
||||
// translation[1] = 0.;
|
||||
// translation[2] = 0.;
|
||||
|
||||
// transform->SetTranslation(translation);
|
||||
|
||||
// // constant for converting degrees to radians
|
||||
// const double dtr = (atan(1.0) * 4.0) / 180.0;
|
||||
// transform->SetRotation(dtr * 0., dtr * 0., dtr * 0.);
|
||||
|
||||
// // The centre of rotation is set by default to the centre of the 3D
|
||||
// // volume but can be offset from this position using a command
|
||||
// // line specified translation [cx,cy,cz]
|
||||
|
||||
// ImageType3D::PointType origin3D = image3DIn->GetOrigin();
|
||||
// const itk::Vector<double, 3> resolution3D = image3DIn->GetSpacing();
|
||||
|
||||
// using ImageRegionType3D = ImageType3D::RegionType;
|
||||
// using SizeType3D = ImageRegionType3D::SizeType;
|
||||
|
||||
// ImageRegionType3D region3D = caster3D->GetOutput()->GetBufferedRegion();
|
||||
// SizeType3D size3D = region3D.GetSize();
|
||||
|
||||
// TransformType::InputPointType isocenter;
|
||||
// if (customized_iso)
|
||||
// {
|
||||
// // // Isocenter location given by the user.
|
||||
// // isocenter[0] = origin3D[0] + resolution3D[0] * cx;
|
||||
// // isocenter[1] = origin3D[1] + resolution3D[1] * cy;
|
||||
// // isocenter[2] = origin3D[2] + resolution3D[2] * cz;
|
||||
// }
|
||||
// else
|
||||
// {
|
||||
// // Set the center of the image as the isocenter.
|
||||
// isocenter[0] = origin3D[0] + resolution3D[0] * static_cast<double>(size3D[0]) / 2.0;
|
||||
// isocenter[1] = origin3D[1] + resolution3D[1] * static_cast<double>(size3D[1]) / 2.0;
|
||||
// isocenter[2] = origin3D[2] + resolution3D[2] * static_cast<double>(size3D[2]) / 2.0;
|
||||
// }
|
||||
|
||||
// transform->SetCenter(isocenter);
|
||||
|
||||
|
||||
// if (true)
|
||||
// {
|
||||
// std::cout << "3D image size: " << size3D[0] << ", " << size3D[1] << ", " << size3D[2] << std::endl
|
||||
// << " resolution: " << resolution3D[0] << ", " << resolution3D[1] << ", " << resolution3D[2] << std::endl
|
||||
// << "Transform: " << transform << std::endl;
|
||||
// }
|
||||
|
||||
|
||||
// // Set the origin of the 2D image
|
||||
// // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
|
||||
|
||||
// // For correct (perspective) projection of the 3D volume, the 2D
|
||||
// // image needs to be placed at a certain distance (the source-to-
|
||||
// // isocenter distance {scd} ) from the focal point, and the normal
|
||||
// // from the imaging plane to the focal point needs to be specified.
|
||||
// //
|
||||
// // By default, the imaging plane normal is set by default to the
|
||||
// // center of the 2D image but may be modified from this using the
|
||||
// // command line parameters [image1centerX, image1centerY,
|
||||
// // image2centerX, image2centerY].
|
||||
|
||||
// double origin2D1[Dimension];
|
||||
// double origin2D2[Dimension];
|
||||
|
||||
// // Note: Two 2D images may have different image sizes and pixel dimensions, although
|
||||
// // scd are the same.
|
||||
|
||||
// const itk::Vector<double, 3> resolution2D1 = rescaler1->GetOutput()->GetSpacing();
|
||||
// const itk::Vector<double, 3> resolution2D2 = rescaler2->GetOutput()->GetSpacing();
|
||||
|
||||
// using ImageRegionType2D = InternalImageType::RegionType;
|
||||
// using SizeType2D = ImageRegionType2D::SizeType;
|
||||
|
||||
// ImageRegionType2D region2D1 = rescaler1->GetOutput()->GetBufferedRegion();
|
||||
// ImageRegionType2D region2D2 = rescaler2->GetOutput()->GetBufferedRegion();
|
||||
// SizeType2D size2D1 = region2D1.GetSize();
|
||||
// SizeType2D size2D2 = region2D2.GetSize();
|
||||
|
||||
// double image1centerX = 0.0;
|
||||
// double image1centerY = 0.0;
|
||||
// double image2centerX = 0.0;
|
||||
// double image2centerY = 0.0;
|
||||
|
||||
// double image1resX = 1.0;
|
||||
// double image1resY = 1.0;
|
||||
// double image2resX = 1.0;
|
||||
// double image2resY = 1.0;
|
||||
|
||||
// if (!customized_2DCX)
|
||||
// { // Central axis positions are not given by the user. Use the image centers
|
||||
// // as the central axis position.
|
||||
// image1centerX = ((double)size2D1[0] - 1.) / 2.;
|
||||
// image1centerY = ((double)size2D1[1] - 1.) / 2.;
|
||||
// image2centerX = ((double)size2D2[0] - 1.) / 2.;
|
||||
// image2centerY = ((double)size2D2[1] - 1.) / 2.;
|
||||
// }
|
||||
|
||||
// // 2D Image 1
|
||||
// origin2D1[0] = -resolution2D1[0] * image1centerX;
|
||||
// origin2D1[1] = -resolution2D1[1] * image1centerY;
|
||||
// origin2D1[2] = -scd;
|
||||
|
||||
// flipFilterLocalizer1->GetOutput()->SetOrigin(origin2D1);
|
||||
// rescaler1->GetOutput()->SetOrigin(origin2D1);
|
||||
|
||||
// // 2D Image 2
|
||||
// origin2D2[0] = -resolution2D2[0] * image2centerX;
|
||||
// origin2D2[1] = -resolution2D2[1] * image2centerY;
|
||||
// origin2D2[2] = -scd;
|
||||
|
||||
// flipFilterLocalizer2->GetOutput()->SetOrigin(origin2D2);
|
||||
// rescaler2->GetOutput()->SetOrigin(origin2D2);
|
||||
|
||||
// registration->SetFixedImageRegion1(rescaler1->GetOutput()->GetBufferedRegion());
|
||||
// registration->SetFixedImageRegion2(rescaler2->GetOutput()->GetBufferedRegion());
|
||||
|
||||
// std::cout<<"BuffRegion ROI 1: "<<rescaler1->GetOutput()->GetBufferedRegion()<<std::endl;
|
||||
// std::cout<<"BuffRegion ROI 2: "<<rescaler2->GetOutput()->GetBufferedRegion()<<std::endl;
|
||||
|
||||
// if (true)
|
||||
// {
|
||||
// std::cout << "2D image 1 size: " << size2D1[0] << ", " << size2D1[1] << ", " << size2D1[2] << std::endl
|
||||
// << " resolution: " << resolution2D1[0] << ", " << resolution2D1[1] << ", " << resolution2D1[2]
|
||||
// << std::endl
|
||||
// << " and position: " << origin2D1[0] << ", " << origin2D1[1] << ", " << origin2D1[2] << std::endl
|
||||
// << "2D image 2 size: " << size2D2[0] << ", " << size2D2[1] << ", " << size2D2[2] << std::endl
|
||||
// << " resolution: " << resolution2D2[0] << ", " << resolution2D2[1] << ", " << resolution2D2[2]
|
||||
// << std::endl
|
||||
// << " and position: " << origin2D2[0] << ", " << origin2D2[1] << ", " << origin2D2[2] << std::endl;
|
||||
// }
|
||||
|
||||
// // Initialize the ray cast interpolator
|
||||
// // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
|
||||
|
||||
// // The ray cast interpolator is used to project the 3D volume. It
|
||||
// // does this by casting rays from the (transformed) focal point to
|
||||
// // each (transformed) pixel coordinate in the 2D image.
|
||||
// //
|
||||
// // In addition a threshold may be specified to ensure that only
|
||||
// // intensities greater than a given value contribute to the
|
||||
// // projected volume. This can be used, for instance, to remove soft
|
||||
// // tissue from projections of CT data and force the registration
|
||||
// // to find a match which aligns bony structures in the images.
|
||||
|
||||
// // 2D Image 1
|
||||
// interpolator1->SetProjectionAngle(dtr * projAngle1);
|
||||
// interpolator1->SetFocalPointToIsocenterDistance(scd);
|
||||
// interpolator1->SetThreshold(threshold);
|
||||
// interpolator1->SetTransform(transform);
|
||||
|
||||
// interpolator1->Initialize();
|
||||
|
||||
// // 2D Image 2
|
||||
// interpolator2->SetProjectionAngle(dtr * projAngle2);
|
||||
// interpolator2->SetFocalPointToIsocenterDistance(scd);
|
||||
// interpolator2->SetThreshold(threshold);
|
||||
// interpolator2->SetTransform(transform);
|
||||
|
||||
// interpolator2->Initialize();
|
||||
|
||||
|
||||
// // Set up the transform and start position
|
||||
// // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
|
||||
|
||||
// // The registration start position is intialised using the
|
||||
// // transformation parameters.
|
||||
|
||||
// registration->SetInitialTransformParameters(transform->GetParameters());
|
||||
|
||||
// // We wish to minimize the negative normalized correlation similarity measure.
|
||||
|
||||
// // optimizer->SetMaximize( true ); // for GradientDifferenceTwoImageToOneImageMetric
|
||||
// optimizer->SetMaximize(false); // for NCC
|
||||
|
||||
// optimizer->SetMaximumIteration(10);
|
||||
// optimizer->SetMaximumLineIteration(4); // for Powell's method
|
||||
// optimizer->SetStepLength(4.0);
|
||||
// optimizer->SetStepTolerance(0.02);
|
||||
// optimizer->SetValueTolerance(0.001);
|
||||
|
||||
// // The optimizer weightings are set such that one degree equates to
|
||||
// // one millimeter.
|
||||
|
||||
// itk::Optimizer::ScalesType weightings(transform->GetNumberOfParameters());
|
||||
|
||||
// weightings[0] = 1. / dtr;
|
||||
// weightings[1] = 1. / dtr;
|
||||
// weightings[2] = 1. / dtr;
|
||||
// weightings[3] = 1.;
|
||||
// weightings[4] = 1.;
|
||||
// weightings[5] = 1.;
|
||||
|
||||
// optimizer->SetScales(weightings);
|
||||
|
||||
// if (true)
|
||||
// {
|
||||
// optimizer->Print(std::cout);
|
||||
// }
|
||||
|
||||
|
||||
// // Create the observers
|
||||
// // ~~~~~~~~~~~~~~~~~~~~
|
||||
|
||||
// CommandIterationUpdate::Pointer observer = CommandIterationUpdate::New();
|
||||
|
||||
// optimizer->AddObserver(itk::IterationEvent(), observer);
|
||||
|
||||
|
||||
// // Start the registration
|
||||
// // ~~~~~~~~~~~~~~~~~~~~~~
|
||||
|
||||
// // Create a timer to record calculation time.
|
||||
// itk::TimeProbesCollectorBase timer;
|
||||
|
||||
// if (true)
|
||||
// {
|
||||
// std::cout << "Starting the registration now..." << std::endl;
|
||||
// }
|
||||
|
||||
// try
|
||||
// {
|
||||
// timer.Start("Registration");
|
||||
// // Start the registration.
|
||||
// registration->StartRegistration();
|
||||
// timer.Stop("Registration");
|
||||
// }
|
||||
// catch (itk::ExceptionObject & err)
|
||||
// {
|
||||
// std::cout << "ExceptionObject caught !" << std::endl;
|
||||
// std::cout << err << std::endl;
|
||||
// return -1;
|
||||
// }
|
||||
|
||||
// using ParametersType = RegistrationType::ParametersType;
|
||||
// ParametersType finalParameters = registration->GetLastTransformParameters();
|
||||
|
||||
// const double RotationAlongX = finalParameters[0] / dtr; // Convert radian to degree
|
||||
// const double RotationAlongY = finalParameters[1] / dtr;
|
||||
// const double RotationAlongZ = finalParameters[2] / dtr;
|
||||
// const double TranslationAlongX = finalParameters[3];
|
||||
// const double TranslationAlongY = finalParameters[4];
|
||||
// const double TranslationAlongZ = finalParameters[5];
|
||||
|
||||
// const int numberOfIterations = optimizer->GetCurrentIteration();
|
||||
|
||||
// const double bestValue = optimizer->GetValue();
|
||||
|
||||
// std::cout << "Result = " << std::endl;
|
||||
// std::cout << " Rotation Along X = " << RotationAlongX << " deg" << std::endl;
|
||||
// std::cout << " Rotation Along Y = " << RotationAlongY << " deg" << std::endl;
|
||||
// std::cout << " Rotation Along Z = " << RotationAlongZ << " deg" << std::endl;
|
||||
// std::cout << " Translation X = " << TranslationAlongX << " mm" << std::endl;
|
||||
// std::cout << " Translation Y = " << TranslationAlongY << " mm" << std::endl;
|
||||
// std::cout << " Translation Z = " << TranslationAlongZ << " mm" << std::endl;
|
||||
// std::cout << " Number Of Iterations = " << numberOfIterations << std::endl;
|
||||
// std::cout << " Metric value = " << bestValue << std::endl;
|
||||
|
||||
|
||||
// // Write out the projection images at the registration position
|
||||
// // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
|
||||
|
||||
// TransformType::Pointer finalTransform = TransformType::New();
|
||||
// // The transform is determined by the parameters and the rotation center.
|
||||
// finalTransform->SetParameters(finalParameters);
|
||||
// finalTransform->SetCenter(isocenter);
|
||||
|
||||
// using ResampleFilterType = itk::ResampleImageFilter<InternalImageType, InternalImageType>;
|
||||
|
||||
// // The ResampleImageFilter is the driving force for the projection image generation.
|
||||
// ResampleFilterType::Pointer resampleFilter1 = ResampleFilterType::New();
|
||||
|
||||
// resampleFilter1->SetInput(caster3D->GetOutput()); // Link the 3D volume.
|
||||
// resampleFilter1->SetDefaultPixelValue(0);
|
||||
|
||||
// // The parameters of interpolator1, such as ProjectionAngle and FocalPointToIsocenterDistance
|
||||
// // have been set before registration. Here we only need to replace the initial
|
||||
// // transform with the final transform.
|
||||
// interpolator1->SetTransform(finalTransform);
|
||||
// interpolator1->Initialize();
|
||||
// resampleFilter1->SetInterpolator(interpolator1);
|
||||
|
||||
// // The output 2D projection image has the same image size, origin, and the pixel spacing as
|
||||
// // those of the input 2D image.
|
||||
// resampleFilter1->SetSize(rescaler1->GetOutput()->GetLargestPossibleRegion().GetSize());
|
||||
// resampleFilter1->SetOutputOrigin(rescaler1->GetOutput()->GetOrigin());
|
||||
// resampleFilter1->SetOutputSpacing(rescaler1->GetOutput()->GetSpacing());
|
||||
|
||||
// // Do the same thing for the output image 2.
|
||||
// ResampleFilterType::Pointer resampleFilter2 = ResampleFilterType::New();
|
||||
// resampleFilter2->SetInput(caster3D->GetOutput());
|
||||
// resampleFilter2->SetDefaultPixelValue(0);
|
||||
|
||||
// // The parameters of interpolator2, such as ProjectionAngle and FocalPointToIsocenterDistance
|
||||
// // have been set before registration. Here we only need to replace the initial
|
||||
// // transform with the final transform.
|
||||
// interpolator2->SetTransform(finalTransform);
|
||||
// interpolator2->Initialize();
|
||||
// resampleFilter2->SetInterpolator(interpolator2);
|
||||
|
||||
// resampleFilter2->SetSize(rescaler2->GetOutput()->GetLargestPossibleRegion().GetSize());
|
||||
// resampleFilter2->SetOutputOrigin(rescaler2->GetOutput()->GetOrigin());
|
||||
// resampleFilter2->SetOutputSpacing(rescaler2->GetOutput()->GetSpacing());
|
||||
|
||||
// /////////////////////////////---DEGUG--START----////////////////////////////////////
|
||||
// if (true)
|
||||
// {
|
||||
// InternalImageType::PointType outputorigin2D1 = rescaler1->GetOutput()->GetOrigin();
|
||||
// std::cout << "Output image 1 origin: " << outputorigin2D1[0] << ", " << outputorigin2D1[1] << ", "
|
||||
// << outputorigin2D1[2] << std::endl;
|
||||
// InternalImageType::PointType outputorigin2D2 = rescaler2->GetOutput()->GetOrigin();
|
||||
// std::cout << "Output image 2 origin: " << outputorigin2D2[0] << ", " << outputorigin2D2[1] << ", "
|
||||
// << outputorigin2D2[2] << std::endl;
|
||||
// }
|
||||
// /////////////////////////////---DEGUG--END----//////////////////////////////////////
|
||||
|
||||
|
||||
// // As explained before, the computed projection is upsided-down.
|
||||
// // Here we use a FilpImageFilter to flip the images in y-direction.
|
||||
// flipFilter1->SetInput(resampleFilter1->GetOutput());
|
||||
// flipFilter2->SetInput(resampleFilter2->GetOutput());
|
||||
|
||||
// // Rescale the intensity of the projection images to 0-255 for output.
|
||||
// //using RescaleFilterType = itk::RescaleIntensityImageFilter<InternalImageType, OutputImageType>;
|
||||
|
||||
// RescaleFilterType::Pointer rescaler1b = RescaleFilterType::New();
|
||||
// rescaler1b->SetOutputMinimum(0);
|
||||
// rescaler1b->SetOutputMaximum(255);
|
||||
// rescaler1b->SetInput(flipFilter1->GetOutput());
|
||||
|
||||
// RescaleFilterType::Pointer rescaler2b = RescaleFilterType::New();
|
||||
// rescaler2b->SetOutputMinimum(0);
|
||||
// rescaler2b->SetOutputMaximum(255);
|
||||
// rescaler2b->SetInput(flipFilter2->GetOutput());
|
||||
|
||||
|
||||
// using WriterType = itk::ImageFileWriter<OutputImageType>;
|
||||
// WriterType::Pointer writer1 = WriterType::New();
|
||||
// WriterType::Pointer writer2 = WriterType::New();
|
||||
|
||||
// writer1->SetFileName("Out1.tif");
|
||||
// writer1->SetInput(0,rescaler1b->GetOutput());
|
||||
|
||||
// try
|
||||
// {
|
||||
// //std::cout << "Writing image: " << fileOutput1 << std::endl;
|
||||
// writer1->Update();
|
||||
// }
|
||||
// catch (itk::ExceptionObject & err)
|
||||
// {
|
||||
// std::cerr << "ERROR: ExceptionObject caught !" << std::endl;
|
||||
// std::cerr << err << std::endl;
|
||||
// }
|
||||
|
||||
// writer2->SetFileName("Out2.tif");
|
||||
// writer2->SetInput(0,rescaler2b->GetOutput());
|
||||
|
||||
// try
|
||||
// {
|
||||
// writer2->Update();
|
||||
// }
|
||||
// catch (itk::ExceptionObject & err)
|
||||
// {
|
||||
// std::cerr << "ERROR: ExceptionObject caught !" << std::endl;
|
||||
// std::cerr << err << std::endl;
|
||||
// }
|
||||
// timer.Report();
|
||||
|
||||
// return 0;
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
//}
|
||||
|
||||
//void itkImageProcessor::RunRegistration()
|
||||
//{
|
||||
|
||||
//}
|
||||
|
||||
|
||||
|
||||
|
||||
@ -1969,7 +1444,7 @@ void itkImageProcessor::SetInitialTranslations(double dX, double dY, double dZ)
|
||||
|
||||
InternalImageType::DirectionType IECtoLPS_Directions;
|
||||
|
||||
|
||||
IECtoLPS_Directions = LPStoIEC_Directions.GetTranspose();
|
||||
|
||||
ImageType3D::PointType LPSTranslations =
|
||||
IECtoLPS_Directions * IECTranslations;
|
||||
@ -1978,9 +1453,9 @@ void itkImageProcessor::SetInitialTranslations(double dX, double dY, double dZ)
|
||||
TZero[1]= LPSTranslations[1];
|
||||
TZero[2]= LPSTranslations[2];
|
||||
|
||||
// TZero[0]= dX;
|
||||
// TZero[1]= dY;
|
||||
// TZero[2]= dZ;
|
||||
|
||||
// std::cout<< "itkImageProcessor::SetInitialTranslations IEC : "<<IECTranslations <<std::endl;
|
||||
// std::cout<< "itkImageProcessor::SetInitialTranslations LPS : "<<TZero[0]<<" "<<TZero[1]<<" "<<TZero[2] <<std::endl;
|
||||
}
|
||||
|
||||
void itkImageProcessor::SetInitialRotations(double dX, double dY, double dZ)
|
||||
|
@ -46,36 +46,36 @@ namespace itk
|
||||
class CommandIterationUpdate : public itk::Command
|
||||
{
|
||||
public:
|
||||
using Self = CommandIterationUpdate;
|
||||
using Superclass = itk::Command;
|
||||
using Pointer = itk::SmartPointer<Self>;
|
||||
itkNewMacro(Self);
|
||||
using Self = CommandIterationUpdate;
|
||||
using Superclass = itk::Command;
|
||||
using Pointer = itk::SmartPointer<Self>;
|
||||
itkNewMacro(Self);
|
||||
|
||||
protected:
|
||||
CommandIterationUpdate() = default;
|
||||
CommandIterationUpdate() = default;
|
||||
|
||||
public:
|
||||
using OptimizerType = itk::PowellOptimizer;
|
||||
using OptimizerPointer = const OptimizerType *;
|
||||
using OptimizerType = itk::PowellOptimizer;
|
||||
using OptimizerPointer = const OptimizerType *;
|
||||
|
||||
void
|
||||
Execute(itk::Object * caller, const itk::EventObject & event) override
|
||||
{
|
||||
Execute((const itk::Object *)caller, event);
|
||||
}
|
||||
|
||||
void
|
||||
Execute(const itk::Object * object, const itk::EventObject & event) override
|
||||
{
|
||||
auto optimizer = dynamic_cast<OptimizerPointer>(object);
|
||||
if (typeid(event) != typeid(itk::IterationEvent))
|
||||
void
|
||||
Execute(itk::Object * caller, const itk::EventObject & event) override
|
||||
{
|
||||
return;
|
||||
Execute((const itk::Object *)caller, event);
|
||||
}
|
||||
|
||||
void
|
||||
Execute(const itk::Object * object, const itk::EventObject & event) override
|
||||
{
|
||||
auto optimizer = dynamic_cast<OptimizerPointer>(object);
|
||||
if (typeid(event) != typeid(itk::IterationEvent))
|
||||
{
|
||||
return;
|
||||
}
|
||||
std::cout << "Iteration: " << optimizer->GetCurrentIteration() << std::endl;
|
||||
std::cout << "Similarity: " << optimizer->GetValue() << std::endl;
|
||||
std::cout << "Position: " << optimizer->GetCurrentPosition() << std::endl;
|
||||
}
|
||||
std::cout << "Iteration: " << optimizer->GetCurrentIteration() << std::endl;
|
||||
std::cout << "Similarity: " << optimizer->GetValue() << std::endl;
|
||||
std::cout << "Position: " << optimizer->GetCurrentPosition() << std::endl;
|
||||
}
|
||||
};
|
||||
|
||||
|
||||
@ -92,215 +92,203 @@ class ITK_EXPORT itkImageProcessor : public Object
|
||||
constexpr static unsigned int Dimension = 3;
|
||||
|
||||
public:
|
||||
/** Standard "Self" typedef. */
|
||||
typedef itkImageProcessor Self;
|
||||
/** Standard "Superclass" typedef. */
|
||||
typedef Object Superclass;
|
||||
/** Smart pointer typedef support */
|
||||
typedef SmartPointer<Self> Pointer;
|
||||
typedef SmartPointer<const Self> ConstPointer;
|
||||
/** Method of creation through the object factory. */
|
||||
itkNewMacro(Self);
|
||||
/** Run-time type information (and related methods). */
|
||||
itkTypeMacro(itkImageProcessor, Object);
|
||||
/** Standard "Self" typedef. */
|
||||
typedef itkImageProcessor Self;
|
||||
/** Standard "Superclass" typedef. */
|
||||
typedef Object Superclass;
|
||||
/** Smart pointer typedef support */
|
||||
typedef SmartPointer<Self> Pointer;
|
||||
typedef SmartPointer<const Self> ConstPointer;
|
||||
/** Method of creation through the object factory. */
|
||||
itkNewMacro(Self);
|
||||
/** Run-time type information (and related methods). */
|
||||
itkTypeMacro(itkImageProcessor, Object);
|
||||
|
||||
/** Input data load methods*/
|
||||
//void load3DVolume(const char *);
|
||||
int load3DSerie(const char* );
|
||||
int load2D(const char *);
|
||||
/** Input data load methods*/
|
||||
//void load3DVolume(const char *);
|
||||
int load3DSerie(const char* );
|
||||
int load2D(const char *);
|
||||
|
||||
/** Projection angles - Gantry angle */
|
||||
void SetProjectionAngle1(double);
|
||||
void SetProjectionAngle2(double);
|
||||
/** Projection angles - Gantry angle */
|
||||
void SetProjectionAngle1(double);
|
||||
void SetProjectionAngle2(double);
|
||||
|
||||
/** Source to axis distance - SAD */
|
||||
void SetSCD(double);
|
||||
/** Source to axis distance - SAD */
|
||||
void SetSCD(double);
|
||||
|
||||
/** Intensity threshold used for image projection,
|
||||
/** Intensity threshold used for image projection,
|
||||
* any voxel below threshold is ignored */
|
||||
void SetIntensityThreshold(double);
|
||||
void SetIntensityThreshold(double);
|
||||
|
||||
/** Custom settings of the projection images */
|
||||
void SetCustom_Isocenter(double,double,double);
|
||||
void SetCustom_2Dres(double,double,double,double);
|
||||
void SetCustom_2Dsize(int,int,int,int);
|
||||
void SetCustom_2Dcenter(double dX1, double dY1, double dX2, double dY2);
|
||||
/** Custom settings of the projection images */
|
||||
void SetCustom_Isocenter(double,double,double);
|
||||
void SetCustom_2Dres(double,double,double,double);
|
||||
void SetCustom_2Dsize(int,int,int,int);
|
||||
void SetCustom_2Dcenter(double dX1, double dY1, double dX2, double dY2);
|
||||
|
||||
/** Set transform parameters for 3D Volume */
|
||||
void SetInitialTranslations(double, double, double);
|
||||
void SetInitialRotations(double, double, double);
|
||||
/** Set transform parameters for 3D Volume */
|
||||
void SetInitialTranslations(double, double, double);
|
||||
void SetInitialRotations(double, double, double);
|
||||
|
||||
/** Get the Localizer intensity window and center for rendering */
|
||||
double* GetImageIntensityWindow(int );
|
||||
/** Get the Localizer intensity window and center for rendering */
|
||||
double* GetImageIntensityWindow(int );
|
||||
|
||||
/** 2D-3D registration routines */
|
||||
// int InitializeRegistration();
|
||||
void InitializeProjector();
|
||||
// void RunRegistration();
|
||||
/** 2D-3D registration routines */
|
||||
// int InitializeRegistration();
|
||||
void InitializeProjector();
|
||||
// void RunRegistration();
|
||||
|
||||
/** Compute Digital Localizers */
|
||||
void GetProjectionImages();
|
||||
/** Compute Digital Localizers */
|
||||
void GetProjectionImages();
|
||||
|
||||
/** Conveniency methods to get VTK images for rendering */
|
||||
vtkImageData* GetProjection1VTK();
|
||||
vtkImageData* GetProjection2VTK();
|
||||
/** Conveniency methods to get VTK images for rendering */
|
||||
vtkImageData* GetProjection1VTK();
|
||||
vtkImageData* GetProjection2VTK();
|
||||
|
||||
vtkImageData* GetLocalizer1VTK();
|
||||
vtkImageData* GetLocalizer2VTK();
|
||||
vtkImageData* GetLocalizer1VTK();
|
||||
vtkImageData* GetLocalizer2VTK();
|
||||
|
||||
/** Debug writers */
|
||||
void WriteProjectionImages();
|
||||
void Write2DImages();
|
||||
/** Debug writers */
|
||||
void WriteProjectionImages();
|
||||
void Write2DImages();
|
||||
|
||||
/** Various pixel types */
|
||||
/** Various pixel types */
|
||||
using InternalPixelType = float;
|
||||
using PixelType2D = short;
|
||||
using PixelType3D = short;
|
||||
using OutputPixelType = unsigned char;
|
||||
|
||||
using ImageType3D = itk::Image<PixelType3D, Dimension>;
|
||||
using InternalImageType = itk::Image<InternalPixelType, Dimension>;
|
||||
typedef enum eProjectionOrientationType
|
||||
{NA = 0, LAT=2, PA=1}
|
||||
tProjOrientationType;
|
||||
/** Enum type for Orientation */
|
||||
using ImageType3D = itk::Image<PixelType3D, Dimension>;
|
||||
using InternalImageType = itk::Image<InternalPixelType, Dimension>;
|
||||
typedef enum eProjectionOrientationType
|
||||
{NA = 0, LAT=2, PA=1}
|
||||
tProjOrientationType;
|
||||
/** Enum type for Orientation */
|
||||
typedef enum eImageOrientationType
|
||||
{ NotDefined = 0, HFS = 1, FFS = 2}
|
||||
tPatOrientation;
|
||||
protected:
|
||||
itkImageProcessor();
|
||||
virtual ~itkImageProcessor();
|
||||
void PrintSelf(std::ostream& os, Indent indent) const;
|
||||
itkImageProcessor();
|
||||
virtual ~itkImageProcessor();
|
||||
void PrintSelf(std::ostream& os, Indent indent) const;
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
private:
|
||||
itkImageProcessor(const Self&); //purposely not implemented
|
||||
void operator=(const Self&); //purposely not implemented
|
||||
itkImageProcessor(const Self&); //purposely not implemented
|
||||
void operator=(const Self&); //purposely not implemented
|
||||
|
||||
|
||||
|
||||
void Finalize();
|
||||
void Finalize();
|
||||
|
||||
|
||||
/** Image types */
|
||||
using ImageType2D = itk::Image<PixelType3D, Dimension>;
|
||||
/** Image types */
|
||||
using ImageType2D = itk::Image<PixelType3D, Dimension>;
|
||||
|
||||
using OutputImageType = itk::Image<OutputPixelType, Dimension>;
|
||||
using OutputImageType = itk::Image<OutputPixelType, Dimension>;
|
||||
|
||||
/** The following lines define each of the components used in the
|
||||
/** The following lines define each of the components used in the
|
||||
registration: The transform, optimizer, metric, interpolator and
|
||||
the registration method itself. */
|
||||
using TransformType = itk::Euler3DTransform<double>;
|
||||
using OptimizerType = itk::PowellOptimizer;
|
||||
using MetricType = itk::NormalizedCorrelationTwoImageToOneImageMetric<InternalImageType, InternalImageType>;
|
||||
using InterpolatorType = itk::gSiddonJacobsRayCastInterpolateImageFunction<InternalImageType, double>;
|
||||
using RegistrationType = itk::TwoProjectionImageRegistrationMethod<InternalImageType, InternalImageType>;
|
||||
using ParametersType = RegistrationType::ParametersType;
|
||||
using TransformType = itk::Euler3DTransform<double>;
|
||||
using OptimizerType = itk::PowellOptimizer;
|
||||
using MetricType = itk::NormalizedCorrelationTwoImageToOneImageMetric<InternalImageType, InternalImageType>;
|
||||
using InterpolatorType = itk::gSiddonJacobsRayCastInterpolateImageFunction<InternalImageType, double>;
|
||||
using RegistrationType = itk::TwoProjectionImageRegistrationMethod<InternalImageType, InternalImageType>;
|
||||
using ParametersType = RegistrationType::ParametersType;
|
||||
|
||||
/** Image reader types */
|
||||
using ImageReaderType2D = itk::ImageFileReader<ImageType2D>;
|
||||
using ImageReaderType3D = itk::ImageFileReader<ImageType3D>;
|
||||
using ImageSeriesReaderType3D = itk::ImageSeriesReader<ImageType3D>;
|
||||
using ImageIOType = itk::GDCMImageIO;
|
||||
using DuplicatorType = itk::ImageDuplicator<InternalImageType>;
|
||||
/** Image reader types */
|
||||
using ImageReaderType2D = itk::ImageFileReader<ImageType2D>;
|
||||
using ImageReaderType3D = itk::ImageFileReader<ImageType3D>;
|
||||
using ImageSeriesReaderType3D = itk::ImageSeriesReader<ImageType3D>;
|
||||
using ImageIOType = itk::GDCMImageIO;
|
||||
using DuplicatorType = itk::ImageDuplicator<InternalImageType>;
|
||||
|
||||
/** widely used filters: flip, cast, rescale, ... */
|
||||
using FlipFilterType = itk::FlipImageFilter<InternalImageType>;
|
||||
using Input2DRescaleFilterType = itk::RescaleIntensityImageFilter<InternalImageType, InternalImageType>;
|
||||
using CastFilterType3D = itk::CastImageFilter<ImageType3D, InternalImageType>;
|
||||
using ChangeInformationFilterType = itk::ChangeInformationImageFilter<InternalImageType>;
|
||||
using ChangeInformationFilterInputType = itk::ChangeInformationImageFilter<ImageType3D>;
|
||||
using ChangeInformationFilterInput2DType = itk::ChangeInformationImageFilter<ImageType2D>;
|
||||
/** widely used filters: flip, cast, rescale, ... */
|
||||
using FlipFilterType = itk::FlipImageFilter<InternalImageType>;
|
||||
using Input2DRescaleFilterType = itk::RescaleIntensityImageFilter<InternalImageType, InternalImageType>;
|
||||
using CastFilterType3D = itk::CastImageFilter<ImageType3D, InternalImageType>;
|
||||
using ChangeInformationFilterType = itk::ChangeInformationImageFilter<InternalImageType>;
|
||||
using ChangeInformationFilterInputType = itk::ChangeInformationImageFilter<ImageType3D>;
|
||||
using ChangeInformationFilterInput2DType = itk::ChangeInformationImageFilter<ImageType2D>;
|
||||
|
||||
/** ITK to VTK filter */
|
||||
using ITKtoVTKFilterType = itk::ImageToVTKImageFilter<OutputImageType>;
|
||||
/** ITK to VTK filter */
|
||||
using ITKtoVTKFilterType = itk::ImageToVTKImageFilter<OutputImageType>;
|
||||
|
||||
|
||||
|
||||
|
||||
MetricType::Pointer metric;
|
||||
TransformType::Pointer transform;
|
||||
OptimizerType::Pointer optimizer;
|
||||
InterpolatorType::Pointer interpolator1;
|
||||
InterpolatorType::Pointer interpolator2;
|
||||
RegistrationType::Pointer registration;
|
||||
TransformType::Pointer transform;
|
||||
|
||||
// ImageReaderType2D::Pointer imageReader2D;
|
||||
ImageReaderType3D::Pointer imageReader3D;
|
||||
ImageSeriesReaderType3D::Pointer imageSeriesReader3D;
|
||||
InterpolatorType::Pointer interpolator1;
|
||||
InterpolatorType::Pointer interpolator2;
|
||||
|
||||
CommandIterationUpdate::Pointer observer;
|
||||
ImageReaderType3D::Pointer imageReader3D;
|
||||
ImageSeriesReaderType3D::Pointer imageSeriesReader3D;
|
||||
|
||||
FlipFilterType::Pointer flipFilter1;
|
||||
FlipFilterType::Pointer flipFilter2;
|
||||
FlipFilterType::Pointer flipFilterLocalizer1;
|
||||
FlipFilterType::Pointer flipFilterLocalizer2;
|
||||
CommandIterationUpdate::Pointer observer;
|
||||
|
||||
CastFilterType3D::Pointer caster3D;
|
||||
CastFilterType3D::Pointer caster2D1;
|
||||
CastFilterType3D::Pointer caster2D2;
|
||||
CastFilterType3D::Pointer caster3D;
|
||||
CastFilterType3D::Pointer caster2D1;
|
||||
CastFilterType3D::Pointer caster2D2;
|
||||
|
||||
|
||||
TransformType::InputPointType isocenter3D;
|
||||
TransformType::InputPointType isocenter3D;
|
||||
|
||||
ITKtoVTKFilterType::Pointer toVTK2D1;
|
||||
ITKtoVTKFilterType::Pointer toVTK2D2;
|
||||
ITKtoVTKFilterType::Pointer toVTKLocalizer1;
|
||||
ITKtoVTKFilterType::Pointer toVTKLocalizer2;
|
||||
ITKtoVTKFilterType::Pointer toVTK2D1;
|
||||
ITKtoVTKFilterType::Pointer toVTK2D2;
|
||||
ITKtoVTKFilterType::Pointer toVTKLocalizer1;
|
||||
ITKtoVTKFilterType::Pointer toVTKLocalizer2;
|
||||
|
||||
ImageType3D::Pointer image3DIn;
|
||||
InternalImageType::Pointer image2D1In;
|
||||
InternalImageType::Pointer image2D2In;
|
||||
ImageType2D::Pointer image2DInGeneric;
|
||||
ImageType3D::Pointer image3DIn;
|
||||
|
||||
InternalImageType::Pointer imageDRT1In;
|
||||
InternalImageType::Pointer imageDRT2In;
|
||||
|
||||
InternalImageType::Pointer TMPImg1;
|
||||
InternalImageType::Pointer TMPImg2;
|
||||
DuplicatorType::Pointer m_LATSourceDupli;
|
||||
DuplicatorType::Pointer m_PASourceDupli;
|
||||
|
||||
DuplicatorType::Pointer m_LATSourceDupli;
|
||||
DuplicatorType::Pointer m_PASourceDupli;
|
||||
|
||||
tPatOrientation
|
||||
tPatOrientation
|
||||
m_3DPatOrientation;
|
||||
InternalImageType::DirectionType
|
||||
InternalImageType::DirectionType
|
||||
LPStoIEC_Directions;
|
||||
|
||||
|
||||
InternalImageType::DirectionType
|
||||
InternalImageType::DirectionType
|
||||
CalcDRTImageDirections(double angle);
|
||||
|
||||
ImageType3D::PointType
|
||||
ImageType3D::PointType
|
||||
CalcDRTImageOrigin(InternalImageType::Pointer m_Image,
|
||||
ImageType3D::PointType m_VolCOOV,
|
||||
double dAngle
|
||||
);
|
||||
double
|
||||
double
|
||||
CalcProjectionAngleLPS(
|
||||
tPatOrientation pOrient,
|
||||
double pAngleIEC);
|
||||
tPatOrientation pOrient,
|
||||
double pAngleIEC);
|
||||
|
||||
double image1res[2], image2res[2];
|
||||
double image1center[2], image2center[2];
|
||||
int image1Size[2], image2Size[2];
|
||||
double isocenter[3];
|
||||
double scd;
|
||||
double threshold;
|
||||
double projAngle1, projAngle2;
|
||||
double projAngle1_IEC, projAngle2_IEC;
|
||||
double image1res[2], image2res[2];
|
||||
double image1center[2], image2center[2];
|
||||
int image1Size[2], image2Size[2];
|
||||
double isocenter[3];
|
||||
double scd;
|
||||
double threshold;
|
||||
double projAngle1, projAngle2;
|
||||
double projAngle1_IEC, projAngle2_IEC;
|
||||
|
||||
double TZero[3], RZero[3];
|
||||
double image3DCOV[3];
|
||||
double TZero[3], RZero[3];
|
||||
double image3DCOV[3];
|
||||
|
||||
bool customized_iso,customized_2DCX,customized_2DRES, customized_2DSIZE;
|
||||
bool image2D1Loaded, image2D2Loaded, image3DLoaded;
|
||||
bool customized_iso,customized_2DCX,customized_2DRES, customized_2DSIZE;
|
||||
bool image2D1Loaded, image2D2Loaded, image3DLoaded;
|
||||
|
||||
double image1IntensityWindow[2], image2IntensityWindow[2];
|
||||
double image1IntensityWindow[2], image2IntensityWindow[2];
|
||||
|
||||
|
||||
InternalImageType::DirectionType Std_DRT2LPS;
|
||||
InternalImageType::DirectionType Std_DRT2LPS;
|
||||
|
||||
};
|
||||
|
||||
|
Reference in New Issue
Block a user