diff --git a/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp b/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp index bc824eb..47c616f 100644 --- a/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp +++ b/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp @@ -25,6 +25,28 @@ gfattori 08.11.2021 #include +/* Some parameters that need to be in settings or read from files */ + +/* Settings */ + +#define CT_IMPORTOFFSET_TXIEC 0 +#define CT_IMPORTOFFSET_TYIEC 0 +#define CT_IMPORTOFFSET_TZIEC 0 +#define CT_IMPORTOFFSET_RXIEC 0 +#define CT_IMPORTOFFSET_RYIEC 0 +#define CT_IMPORTOFFSET_RZIEC 0 + +#define CT_DEFAULT_TABLEHEIGHT 175 + +/* From Files */ +#define RT_iso_x 0 +#define RT_iso_y 0 +#define RT_iso_z 0 +#define RT_LAT 0 +#define RT_LNG 0 +#define RT_VRT 0 + + namespace itk { @@ -32,9 +54,9 @@ namespace itk * but we keep this for ourselves... */ static double Standard_DRT2LPS [9] = { - 1,0,0, - 0,0,-1, - 0,1,0 }; + 1,0,0, + 0,0,-1, + 0,1,0 }; @@ -46,7 +68,7 @@ static double PAElementsIEC[9] = { static double LATElementsIEC[9] = { 0, 0, -1, 0, -1, 0, - -1, 0, 0}; + -1, 0, 0}; static double HFS2IEC[9] = { 1, 0, 0, @@ -84,7 +106,7 @@ const char *gGetStringValueFromTag(const gdcm::Tag& t, const gdcm::DataSet& ds) itkImageProcessor::itkImageProcessor() { -// std::cout<<"itkImageProcessor::itkImageProcessor contructor"<c_str(); - seriesItr++; + std::string seriesIdentifier; + seriesIdentifier = seriesItr->c_str(); + seriesItr++; - std::cout << "\nReading: "; - std::cout << seriesIdentifier << std::endl; - using FileNamesContainer = std::vector; - FileNamesContainer fileNames = nameGenerator->GetFileNames(seriesIdentifier); + std::cout << "\nReading: "; + std::cout << seriesIdentifier << std::endl; + using FileNamesContainer = std::vector; + FileNamesContainer fileNames = nameGenerator->GetFileNames(seriesIdentifier); - using ImageIOType = itk::GDCMImageIO; - ImageIOType::Pointer dicomIO = ImageIOType::New(); - imageSeriesReader3D->SetImageIO(dicomIO); - imageSeriesReader3D->SetFileNames(fileNames); - imageSeriesReader3D->SetNumberOfWorkUnits(8); - imageSeriesReader3D->ForceOrthogonalDirectionOff(); // properly read CTs with gantry tilt + using ImageIOType = itk::GDCMImageIO; + ImageIOType::Pointer dicomIO = ImageIOType::New(); + imageSeriesReader3D->SetImageIO(dicomIO); + imageSeriesReader3D->SetFileNames(fileNames); + imageSeriesReader3D->SetNumberOfWorkUnits(8); + imageSeriesReader3D->ForceOrthogonalDirectionOff(); // properly read CTs with gantry tilt - imageSeriesReader3D->Update(); + imageSeriesReader3D->Update(); - eImageOrientationType - currImgOrient; - /* check patient orientation */ - if(fileNames.size()>0){ +// eImageOrientationType +// currImgOrient; + /* check patient orientation */ + if(fileNames.size()>0){ - gdcm::Reader R; - R.SetFileName(fileNames.at(0).c_str()); - if(!R.Read()) - { - cerr<<"ERROR: cannot read file: "<GetOutput(); - /* before changing origin for projection, let's save the center of volume */ - using ImageRegionType3D = ImageType3D::RegionType; - using SizeType3D = ImageRegionType3D::SizeType; - using ImageDirectionType3D = ImageType3D::DirectionType; + image3DIn = imageSeriesReader3D->GetOutput(); - ImageRegionType3D region3D = image3DIn->GetBufferedRegion(); - SizeType3D size3D = region3D.GetSize(); - ImageType3D::PointType origin3D = image3DIn->GetOrigin(); - const itk::Vector resolution3D = image3DIn->GetSpacing(); - ImageDirectionType3D imagDirection = image3DIn->GetDirection(); + + /* TODO: here is the right place to change origin when required + * after this point, reference info are saved */ + + /* End of origin change */ + + + /* before changing origin for projection, let's save the center of volume */ + using ImageRegionType3D = ImageType3D::RegionType; + using SizeType3D = ImageRegionType3D::SizeType; + using ImageDirectionType3D = ImageType3D::DirectionType; + + ImageRegionType3D region3D = image3DIn->GetBufferedRegion(); + SizeType3D size3D = region3D.GetSize(); +// ImageType3D::PointType origin3D = image3DIn->GetOrigin(); + const itk::Vector resolution3D = image3DIn->GetSpacing(); + ImageDirectionType3D imagDirection = image3DIn->GetDirection(); + + /* 3D volume origin in patient coordinates */ + m_3DOriginLPS = image3DIn->GetOrigin(); + + + /* calculate image size in 3D Space by finding the last voxel position */ + using VectorType = itk::Vector; + 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 = + m_3DOriginLPS + (imagDirection * Dest); + + ImageType3D::PointType image3DCOVLPS; + image3DCOVLPS[0] = (m_3DOriginLPS[0]+LastVoxelPosition[0])/2; + image3DCOVLPS[1] = (m_3DOriginLPS[1]+LastVoxelPosition[1])/2; + image3DCOVLPS[2] = (m_3DOriginLPS[2]+LastVoxelPosition[2])/2; + + /* Calculate center of projection in patient coordinates */ + + /* center of projection in an IEC reference at Patient Origin */ + ImageType3D::PointType ProjectionCenterFixedAxes; + ProjectionCenterFixedAxes[0] = 0.; + ProjectionCenterFixedAxes[1] = 0.; + ProjectionCenterFixedAxes[2] = CT_DEFAULT_TABLEHEIGHT; + + + InternalImageType::DirectionType IECtoLPS_Directions; + IECtoLPS_Directions = LPStoIEC_Directions.GetTranspose(); + + m_3DProjectionOriginLPS = + IECtoLPS_Directions * ProjectionCenterFixedAxes ; + /* in longitudinal direction (Sup-Inf), we put it in the + * middle of the volume */ + m_3DProjectionOriginLPS[2] = image3DCOVLPS[2]; + + /* calculate also the same center with zero origin + * This is required to bring back the DRT into LPS */ + m_3DProjectionOriginLPSZero = + m_3DProjectionOriginLPS - m_3DOriginLPS; - /* calculate image size in 3D Space by fiding the last voxel position */ - using VectorType = itk::Vector; - VectorType Dest; - Dest[0]=(size3D[0]-1) * resolution3D[0]; - Dest[1]=(size3D[1]-1) * resolution3D[1]; - Dest[2]=(size3D[2]-1) * resolution3D[2]; + if(true){ + std::cout<<" -------- 3D IMAGE --------"<SetOrigin(image3DOrigin); + // To simply Siddon-Jacob's fast ray-tracing algorithm, we force the origin of the CT image + // to be (0,0,0). Because we align the CT isocenter with the central axis, the projection + // geometry is fully defined. The origin of the CT image becomes irrelavent. + ImageType3D::PointType image3DOrigin; + image3DOrigin[0] = 0.0; + image3DOrigin[1] = 0.0; + image3DOrigin[2] = 0.0; + image3DIn->SetOrigin(image3DOrigin); - caster3D->SetInput(image3DIn); - caster3D->Update(); + caster3D->SetInput(image3DIn); + caster3D->Update(); } @@ -405,80 +470,68 @@ int itkImageProcessor::load3DSerie(const char * pcDirName) } double - itkImageProcessor::CalcProjectionAngleLPS( +itkImageProcessor::CalcProjectionAngleLPS( tPatOrientation pOrient, double pAngleIEC){ -// std::cout<< "CalcProjectèionAngleLPS :: pAngleIEC " << pAngleIEC <SetInput(imageReader2D->GetOutput()); - caster2DInput->Update(); + /* cast type to internal processing type */ + CastFilterType3D::Pointer caster2DInput; + caster2DInput = CastFilterType3D::New(); + caster2DInput->SetInput(imageReader2D->GetOutput()); + caster2DInput->Update(); - InternalImageType::Pointer MyLocalizerImage = - caster2DInput->GetOutput(); + InternalImageType::Pointer MyLocalizerImage = + caster2DInput->GetOutput(); - ImageType3D::PointType image3DOrigin; - image3DOrigin[0]=image3DCOV[0]; - image3DOrigin[1]=image3DCOV[1]; - image3DOrigin[2]=image3DCOV[2]; +// ChangeInformationFilterType::Pointer filter = ChangeInformationFilterType::New(); +// filter->SetInput(MyLocalizerImage); - ImageType3D::PointType image2Dcenter_3DCOV; - InternalImageType::DirectionType image2D1InDirectionInv; - image2D1InDirectionInv.GetVnlMatrix() = MyLocalizerImage->GetDirection().GetTranspose(); +// try +// { +// filter->UpdateOutputInformation(); +// } +// catch (itk::ExceptionObject & error) +// { +// std::cerr << "Error: " << error << std::endl; +// return EXIT_FAILURE; +// } - image2Dcenter_3DCOV = - - ( - image2D1InDirectionInv * (image3DOrigin - MyLocalizerImage->GetOrigin()) - ); - - /* - We have to set the z coordinate to the focal distance, it would be 0 otherwise... - */ - image2Dcenter_3DCOV[2] = -scd; - - ChangeInformationFilterType::Pointer filter = ChangeInformationFilterType::New(); - filter->SetInput(MyLocalizerImage); - - const InternalImageType::DirectionType newDirection = - MyLocalizerImage->GetDirection() * image2D1InDirectionInv; - - try - { - filter->UpdateOutputInformation(); - } - catch (itk::ExceptionObject & error) - { - std::cerr << "Error: " << error << std::endl; - return EXIT_FAILURE; - } - - filter->Update(); +// filter->Update(); - double* m_ImageIntensity; - FlipFilterType::Pointer m_FlipFilter; - DuplicatorType::Pointer m_Duplicator; - bool* m_ImLoaded; + double* m_ImageIntensity; + FlipFilterType::Pointer m_FlipFilter; + DuplicatorType::Pointer m_Duplicator; + bool* m_ImLoaded; - switch (currProjOrientation) { - case eProjectionOrientationType::PA: + switch (currProjOrientation) { + case eProjectionOrientationType::PA: - m_ImageIntensity = image1IntensityWindow; - m_Duplicator = m_PASourceDupli; - m_ImLoaded = &image2D1Loaded; + m_ImageIntensity = image1IntensityWindow; + m_Duplicator = m_PASourceDupli; +// m_ImLoaded = &image2D1Loaded; - break; + break; - case eProjectionOrientationType::LAT: + case eProjectionOrientationType::LAT: - m_ImageIntensity = image2IntensityWindow; - m_Duplicator = m_LATSourceDupli; - m_ImLoaded = &image2D2Loaded; + m_ImageIntensity = image2IntensityWindow; + m_Duplicator = m_LATSourceDupli; +// m_ImLoaded = &image2D2Loaded; - break; + break; - default: - return -1; - break; - } + default: + return -1; + break; + } - memcpy(m_ImageIntensity,dIntensityWindow, 2*sizeof(double)); + memcpy(m_ImageIntensity,dIntensityWindow, 2*sizeof(double)); - m_Duplicator->SetInputImage(filter->GetOutput()); - m_Duplicator->Update(); - -// m_FlipFilter->SetInput(m_Image); + m_Duplicator->SetInputImage(caster2DInput->GetOutput() );//filter->GetOutput()); + m_Duplicator->Update(); - return - (int) currProjOrientation; + return + (int) currProjOrientation; } @@ -768,39 +797,15 @@ void itkImageProcessor::InitializeProjector() const double dtr = (atan(1.0) * 4.0) / 180.0; transform->SetRotation(dtr * RZero[0], dtr * RZero[1], dtr * RZero[2]); - // 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] + isocenter3D = m_3DProjectionOriginLPSZero; - ImageType3D::PointType origin3D = caster3D->GetOutput()->GetOrigin(); - const itk::Vector resolution3D = caster3D->GetOutput()->GetSpacing(); - using ImageRegionType3D = ImageType3D::RegionType; - using SizeType3D = ImageRegionType3D::SizeType; + // std::cout<<" -------- PROJECTOR GEOMETRY --------"<(size3D[0])<<" "<(size3D[0]) / 2.0<GetOutput()->GetBufferedRegion(); - SizeType3D size3D = region3D.GetSize(); - - if (customized_iso) - { - // Isocenter location given by the user. - 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"<(size3D[0]) / 2.0; - isocenter3D[1] = origin3D[1] + resolution3D[1] * static_cast(size3D[1]) / 2.0; - isocenter3D[2] = origin3D[2] + resolution3D[2] * static_cast(size3D[2]) / 2.0; - } - -// std::cout<<" -------- PROJECTOR GEOMETRY --------"<SetCenter(isocenter3D); // 2D Image 1 @@ -823,23 +828,24 @@ void itkImageProcessor::InitializeProjector() #include + itkImageProcessor::InternalImageType::DirectionType - itkImageProcessor::CalcDRTImageDirections(double angle){ +itkImageProcessor::CalcDRTImageDirections(double angle){ InternalImageType::DirectionType DRTDirs; - DRTDirs.SetIdentity(); - DRTDirs. GetVnlMatrix()[0][0] = cos (angle*itk::Math::pi/180); - DRTDirs. GetVnlMatrix()[0][1] = -sin (angle*itk::Math::pi/180); - DRTDirs. GetVnlMatrix()[1][0] = sin (angle*itk::Math::pi/180); - DRTDirs. GetVnlMatrix()[1][1] = cos (angle*itk::Math::pi/180); + DRTDirs.SetIdentity(); + DRTDirs. GetVnlMatrix()[0][0] = cos (angle*itk::Math::pi/180); + DRTDirs. GetVnlMatrix()[0][1] = -sin (angle*itk::Math::pi/180); + DRTDirs. GetVnlMatrix()[1][0] = sin (angle*itk::Math::pi/180); + DRTDirs. GetVnlMatrix()[1][1] = cos (angle*itk::Math::pi/180); - DRTDirs *= Std_DRT2LPS; + DRTDirs *= Std_DRT2LPS; - return - DRTDirs; + return + DRTDirs; } @@ -850,47 +856,47 @@ itkImageProcessor::InternalImageType::DirectionType */ itkImageProcessor::ImageType3D::PointType - itkImageProcessor::CalcDRTImageOrigin(itkImageProcessor::InternalImageType::Pointer m_Image, - itkImageProcessor::ImageType3D::PointType m_VolCOOV, - double dAngle - ){ +itkImageProcessor::CalcDRTImageOrigin(itkImageProcessor::InternalImageType::Pointer m_Image, + itkImageProcessor::ImageType3D::PointType m_VolCOOV, + double dAngle + ){ - using ImageRegionType3D = itkImageProcessor::ImageType3D::RegionType; - using SizeType3D = ImageRegionType3D::SizeType; - using ImageDirectionType3D = itkImageProcessor::ImageType3D::DirectionType; - using VectorType = itk::Vector; + using ImageRegionType3D = itkImageProcessor::ImageType3D::RegionType; + using SizeType3D = ImageRegionType3D::SizeType; + using ImageDirectionType3D = itkImageProcessor::ImageType3D::DirectionType; + using VectorType = itk::Vector; - ImageRegionType3D region3D; - SizeType3D size3D; - itk::Vector resolution3D; + ImageRegionType3D region3D; + SizeType3D size3D; + itk::Vector resolution3D; - region3D = m_Image ->GetBufferedRegion(); - size3D = region3D.GetSize(); - resolution3D = m_Image ->GetSpacing(); + region3D = m_Image ->GetBufferedRegion(); + size3D = region3D.GetSize(); + resolution3D = m_Image ->GetSpacing(); - VectorType Dest; - Dest[0]=(size3D[0]-1) * resolution3D[0]; - Dest[1]=(size3D[1]-1) * resolution3D[1]; - Dest[2]=(size3D[2]-1) * resolution3D[2]; + VectorType Dest; + Dest[0]=(size3D[0]-1) * resolution3D[0]; + Dest[1]=(size3D[1]-1) * resolution3D[1]; + Dest[2]=(size3D[2]-1) * resolution3D[2]; - itkImageProcessor::ImageType3D::PointType LastVoxelPosition = - m_Image ->GetOrigin() + - (m_Image ->GetDirection() * Dest); + itkImageProcessor::ImageType3D::PointType LastVoxelPosition = + m_Image ->GetOrigin() + + (m_Image ->GetDirection() * Dest); - itkImageProcessor::ImageType3D::PointType DRTCOV; - DRTCOV[0] = (( m_Image ->GetOrigin())[0]+LastVoxelPosition[0])/2; - DRTCOV[1] = (( m_Image ->GetOrigin())[1]+LastVoxelPosition[1])/2; - DRTCOV[2] = (( m_Image ->GetOrigin())[2]+LastVoxelPosition[2])/2; + itkImageProcessor::ImageType3D::PointType DRTCOV; + DRTCOV[0] = (( m_Image ->GetOrigin())[0]+LastVoxelPosition[0])/2; + DRTCOV[1] = (( m_Image ->GetOrigin())[1]+LastVoxelPosition[1])/2; + DRTCOV[2] = (( m_Image ->GetOrigin())[2]+LastVoxelPosition[2])/2; - itkImageProcessor::InternalImageType::DirectionType DRT2LPS - = this->CalcDRTImageDirections(dAngle); + itkImageProcessor::InternalImageType::DirectionType DRT2LPS + = this->CalcDRTImageDirections(dAngle); - ImageType3D::PointType NewOrigin = - m_VolCOOV + DRT2LPS * (m_Image ->GetOrigin() - DRTCOV); + ImageType3D::PointType NewOrigin = + m_VolCOOV + DRT2LPS * (m_Image ->GetOrigin() - DRTCOV); - return - NewOrigin; + return + NewOrigin; } @@ -935,9 +941,9 @@ void itkImageProcessor::GetProjectionImages(){ if(false){ // if(image2D1Loaded){ -// resampleFilter1->SetSize(m_Input2DImage_PA ->GetLargestPossibleRegion().GetSize()); -// resampleFilter1->SetOutputOrigin(m_Input2DImage_PA ->GetOrigin()); -// resampleFilter1->SetOutputSpacing(m_Input2DImage_PA ->GetSpacing()); + // resampleFilter1->SetSize(m_Input2DImage_PA ->GetLargestPossibleRegion().GetSize()); + // resampleFilter1->SetOutputOrigin(m_Input2DImage_PA ->GetOrigin()); + // resampleFilter1->SetOutputSpacing(m_Input2DImage_PA ->GetSpacing()); // std::cout<<"Size "<GetOutput()->GetLargestPossibleRegion().GetSize()<GetOutput()->GetOrigin()<SetSize(m_Input2DImage_LAT->GetLargestPossibleRegion().GetSize()); -// resampleFilter2->SetOutputOrigin(m_Input2DImage_LAT ->GetOrigin()); -// resampleFilter2->SetOutputSpacing(m_Input2DImage_LAT ->GetSpacing()); + // resampleFilter2->SetSize(m_Input2DImage_LAT->GetLargestPossibleRegion().GetSize()); + // resampleFilter2->SetOutputOrigin(m_Input2DImage_LAT ->GetOrigin()); + // resampleFilter2->SetOutputSpacing(m_Input2DImage_LAT ->GetSpacing()); } else { @@ -1026,19 +1032,19 @@ void itkImageProcessor::GetProjectionImages(){ /* First of all we need the center of the Projection images in its reference frame */ - ImageType3D::PointType image3DOrigin; - image3DOrigin[0]=image3DCOV[0]; - image3DOrigin[1]=image3DCOV[1]; - image3DOrigin[2]=image3DCOV[2]; + ImageType3D::PointType image3DOrigin = m_3DProjectionOriginLPS; +// image3DOrigin[0]=image3DCOV[0]; +// image3DOrigin[1]=image3DCOV[1]; +// image3DOrigin[2]=image3DCOV[2]; - resampleFilter1->Update(); + resampleFilter1->Update(); ImageType3D::PointType NewOrigin = this->CalcDRTImageOrigin(resampleFilter1->GetOutput(), - image3DOrigin, + m_3DProjectionOriginLPS, projAngle1); InternalImageType::DirectionType DRT2LPS1 = - this->CalcDRTImageDirections(projAngle1); + this->CalcDRTImageDirections(projAngle1); filter1->SetInput( resampleFilter1->GetOutput() ); filter1->SetOutputDirection(DRT2LPS1 );//IECtoLPS_Directions); @@ -1049,20 +1055,20 @@ void itkImageProcessor::GetProjectionImages(){ - InternalImageType::DirectionType LATtoHFS = - this->CalcDRTImageDirections(projAngle2); +// InternalImageType::DirectionType LATtoHFS = +// this->CalcDRTImageDirections(projAngle2); /* Again.. */ resampleFilter2->Update(); NewOrigin = - this->CalcDRTImageOrigin(resampleFilter2->GetOutput(), - image3DOrigin, - projAngle2); + this->CalcDRTImageOrigin(resampleFilter2->GetOutput(), + m_3DProjectionOriginLPS, + projAngle2); - DRT2LPS1 = - this->CalcDRTImageDirections(projAngle2); + DRT2LPS1 = + this->CalcDRTImageDirections(projAngle2); @@ -1239,19 +1245,19 @@ vtkImageData* itkImageProcessor::GetLocalizer2VTK() ImageType3D::PointType LastVoxelPosition = origin3D + (imagDirection * Dest); -// std::cout<<" -------- Localizer 2 ITK --------"<GetOutput()->GetBounds(); -// std::cout<< "-------- Localizer 2 VTK --------" <GetOutput()->GetBounds(); + // std::cout<< "-------- Localizer 2 VTK --------" <SetOutputMinimum(0); rescaler1->SetOutputMaximum(255); rescaler1->SetInput( imageDRT1In ); -// flipFilter1->GetOutput()); + // flipFilter1->GetOutput()); rescaler1->Update(); @@ -1278,41 +1284,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 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 resolution3D = rescaler1->GetOutput()->GetSpacing(); + // ImageDirectionType3D imagDirection = rescaler1->GetOutput()->GetDirection(); -// /* calculate image size in 3D Space */ -// using VectorType = itk::Vector; -// 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; + // 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 --------"<GetOutput()->GetBounds(); + // double* dBounds = toVTK2D1->GetOutput()->GetBounds(); -// std::cout<< "-------- Proj 1 --------" <GetOutput(); @@ -1327,50 +1333,50 @@ vtkImageData* itkImageProcessor::GetProjection2VTK() rescaler2->SetOutputMinimum(0); rescaler2->SetOutputMaximum(255); rescaler2->SetInput( imageDRT2In ); -// flipFilter2->GetOutput()); + // 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 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 resolution3D = rescaler2->GetOutput()->GetSpacing(); + // ImageDirectionType3D imagDirection = rescaler2->GetOutput()->GetDirection(); /* calculate image size in 3D Space */ -// using VectorType = itk::Vector; -// 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; + // 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 --------"<SetInput(rescaler2->GetOutput()); toVTK2D2->Update(); -// double* dBounds = toVTK2D2->GetOutput()->GetBounds(); + // double* dBounds = toVTK2D2->GetOutput()->GetBounds(); -// std::cout<< "-------- Proj 2 --------" <GetOutput(); @@ -1386,13 +1392,13 @@ void itkImageProcessor::WriteProjectionImages() rescaler1->SetOutputMinimum(0); rescaler1->SetOutputMaximum(255); rescaler1->SetInput( imageDRT1In ); -// flipFilter1->GetOutput()); + // flipFilter1->GetOutput()); RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New(); rescaler2->SetOutputMinimum(0); rescaler2->SetOutputMaximum(255); rescaler2->SetInput( imageDRT2In ); -// flipFilter2->GetOutput()); + // flipFilter2->GetOutput()); using WriterType = itk::ImageFileWriter; @@ -1454,8 +1460,8 @@ void itkImageProcessor::SetInitialTranslations(double dX, double dY, double dZ) TZero[2]= LPSTranslations[2]; -// std::cout<< "itkImageProcessor::SetInitialTranslations IEC : "<