diff --git a/reg23Topograms/itkDTRrecon/DRTMetaInformation.cpp b/reg23Topograms/itkDTRrecon/DRTMetaInformation.cpp index 40d6ff2..c8ab0be 100644 --- a/reg23Topograms/itkDTRrecon/DRTMetaInformation.cpp +++ b/reg23Topograms/itkDTRrecon/DRTMetaInformation.cpp @@ -314,6 +314,39 @@ RTGeometryMetaInformation +TransformMetaInformation :: +TransformMetaInformation (){ + + m_Translations.Fill(0.); + + m_Rotations.Fill(0.); + + m_ReferenceTransform.Fill(0.); + + m_CurrentTransform.Fill(0.); + +} + +void +TransformMetaInformation +::PrintSelf(std::ostream& os, itk::Indent indent) const +{ + Superclass::PrintSelf(os, indent); +} + + +TransformMetaInformation +::~TransformMetaInformation () +{ + +} + + + + + + + } diff --git a/reg23Topograms/itkDTRrecon/DRTMetaInformation.h b/reg23Topograms/itkDTRrecon/DRTMetaInformation.h index d997ee0..c98ea86 100644 --- a/reg23Topograms/itkDTRrecon/DRTMetaInformation.h +++ b/reg23Topograms/itkDTRrecon/DRTMetaInformation.h @@ -217,7 +217,6 @@ public: itkSetMacro(ImageDirections,DirectionType); itkGetMacro(ImageDirections,DirectionType); - //itkSetMacro(LPS2IECDirections,DirectionType); itkGetMacro(LPS2IECDirections,DirectionType); PointType GetCOV(); @@ -328,7 +327,10 @@ protected: m_DRT1Spacing, m_DRT2Spacing; + PointType + /* Transform between IEC Support and + * IEC Scanner frame of reference */ m_IECS2IECScannerT, m_IECS2IECScannerR; @@ -353,7 +355,6 @@ private: - class RTGeometryMetaInformation : public itk::Object{ @@ -384,8 +385,12 @@ public: protected: + PointType + /* RT Plan isocenter in LPS coordinates */ m_IsocenterLPS, + /* RT Plan couch SETUP offset corresponding to + * the RT Plan isocenter */ m_IsocenterIECS; /** Default Constructor **/ @@ -401,6 +406,53 @@ private: void operator=(const Self&); }; + +class TransformMetaInformation : + public itk::Object{ + +public: + /** standard typedefs **/ + typedef TransformMetaInformation Self; + typedef itk::Object Superclass; + typedef itk::SmartPointer Pointer; + typedef itk::Point PointType; + typedef itk::Matrix TransformMatrixType; + + /** Method for creation through the object factory. */ + itkNewMacro(Self); + + /** Run-time type information (and related methods). */ + itkTypeMacro(TransformMetaInformation, itk::Object); + + /** object information streaming **/ + void PrintSelf(std::ostream& os, itk::Indent indent) const; + + +protected: + + + PointType + m_Translations, + m_Rotations; + + TransformMatrixType + m_ReferenceTransform, + m_CurrentTransform; + + + /** Default Constructor **/ + TransformMetaInformation (); + /** Default Destructor **/ + virtual ~TransformMetaInformation (); + +private: + /** purposely not implemented **/ + TransformMetaInformation (const Self&); + + /** purposely not implemented **/ + void operator=(const Self&); +}; + } #endif diff --git a/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp b/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp index 15e3ec0..919aa82 100644 --- a/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp +++ b/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp @@ -25,12 +25,6 @@ gfattori 08.11.2021 #include -/* Some parameters that need to be in settings or read from files */ - -/* Settings */ - -#define CT_DEFAULT_TABLEHEIGHT 175 - namespace itk @@ -44,8 +38,6 @@ static double Standard_DRT2LPS [9] = { 0,0,-1, 0,1,0 }; - - static double PAElementsIEC[9] = { 1, 0, 0, 0, -1, 0, @@ -97,31 +89,14 @@ itkImageProcessor::itkImageProcessor() interpolator1 = InterpolatorType::New(); interpolator2 = InterpolatorType::New(); - -// image1res[0]=image2res[1]= 2.; -// image1Size[0] = image1Size[1] = 512; -// image2Size[0] = image2Size[1] = 512; - - TZero[0]=TZero[1]=TZero[2]=0.; RZero[0]=RZero[1]=RZero[2]=0.; - customized_2DRES = false; - customized_2DSIZE = false; - customized_ImportTransform = false; - customized_ProjectionCenter = false; - customized_RTCouchSetup = false; - customized_RTIsocenter = false; - - image2D1Loaded= false; - image2D2Loaded= false; - image3DLoaded= false; toVTK2D1 = ITKtoVTKFilterType::New(); toVTK2D2 = ITKtoVTKFilterType::New(); toVTKLocalizer1 = ITKtoVTKFilterType::New(); toVTKLocalizer2 = ITKtoVTKFilterType::New(); -// m_3DPatOrientation = eImageOrientationType::NotDefined; m_LATSourceDupli = DuplicatorType::New(); m_PASourceDupli = DuplicatorType::New(); @@ -137,12 +112,6 @@ itkImageProcessor::itkImageProcessor() m_3DInputChangeInformationToZero = ChangeInformationFilterType::New(); -// IEC2DCMMapT.Fill(0.); -// IEC2DCMMapR.Fill(0.); -// rtIsocenterLPS.Fill(0.); -// rtCouchOffset.Fill(0.); - //ProjectionCenterFixedAxes.Fill(0.); - m_Projection1VTKTransform = vtkMatrix4x4::New(); m_Projection1VTKTransform->Identity(); @@ -162,6 +131,11 @@ itkImageProcessor::itkImageProcessor() m_RTMetaInfo = NULL; + m_DRTGeometryMetaInfo = NULL; + + m_TransformMetaInfo = NULL; + + /* Initialiser the projection geoemtry with defaults */ m_DRTGeometryMetaInfo = DRTProjectionGeometryImageMetaInformation::New(); @@ -195,9 +169,6 @@ itkImageProcessor::itkImageProcessor() m_DRTGeometryMetaInfo->SetIECS2IECScannerT(Point3D); m_DRTGeometryMetaInfo->SetIECS2IECScannerR(Point3D); - - - } itkImageProcessor::~itkImageProcessor() @@ -213,7 +184,6 @@ void itkImageProcessor::PrintSelf(std::ostream& os, Indent indent) const void itkImageProcessor::SetIntensityThreshold(double dT) { m_DRTGeometryMetaInfo->SetIntensityThreshold(dT); -// d_intensityThreshold = dT; } void itkImageProcessor::SetSCD(double dDist) @@ -233,13 +203,6 @@ void itkImageProcessor::SetCustom_ImportTransform(double dTx, double dRy, double dRz) { -// IEC2DCMMapT[0] = dTx; -// IEC2DCMMapT[1] = dTy; -// IEC2DCMMapT[2] = dTz; -// IEC2DCMMapR[0] = dRx; -// IEC2DCMMapR[1] = dRy; -// IEC2DCMMapR[2] = dRz; - customized_ImportTransform = true; ImageType3D::PointType Punto; @@ -257,29 +220,9 @@ void itkImageProcessor::SetCustom_ImportTransform(double dTx, } -//void itkImageProcessor::SetCustom_RTCouchSetup_IEC(double dLAT, -// double dVRT, -// double dLNG) -//{ -// rtCouchOffset[0] = dLAT; -// rtCouchOffset[1] = dVRT; -// rtCouchOffset[2] = dLNG; -// customized_RTCouchSetup = true; -//} - -//void itkImageProcessor::SetCustom_RTIsocenter_LPS(double dX, -// double dY, -// double dZ) -//{ -// rtIsocenterLPS[0] = dX; -// rtIsocenterLPS[1] = dY; -// rtIsocenterLPS[2] = dZ; -// customized_RTIsocenter = true; -//} - void itkImageProcessor::SetCustom_ProjectionCenterFixedAxes_IEC(double dX, - double dY, - double dZ) + double dY, + double dZ) { typedef itk::Point PointType; @@ -291,16 +234,11 @@ void itkImageProcessor::SetCustom_ProjectionCenterFixedAxes_IEC(double dX, m_DRTGeometryMetaInfo->SetProjectionCenter(m_point); - customized_ProjectionCenter = true; } void itkImageProcessor::SetCustom_2Dres(double nX1,double nY1,double nX2,double nY2) { -// image1res[0] = nX1; -// image1res[1] = nY1; -// image2res[0] = nX2; -// image2res[1] = nY2; - customized_2DRES = true; + if(m_DRTGeometryMetaInfo == NULL) { // todo @@ -309,23 +247,18 @@ void itkImageProcessor::SetCustom_2Dres(double nX1,double nY1,double nX2,double ImageType3D::SpacingType Spacing; Spacing [0] = nX1; Spacing [1] = nY1; - Spacing [1] = 1.; + Spacing [2] = 1.; m_DRTGeometryMetaInfo->SetDRT1Spacing(Spacing); Spacing [0] = nX2; Spacing [1] = nY2; - Spacing [1] = 1.; + Spacing [2] = 1.; m_DRTGeometryMetaInfo->SetDRT2Spacing(Spacing); //TODO UPDATE TO FOLLOW } void itkImageProcessor::SetCustom_2Dsize(int nX1, int nY1,int nX2,int nY2) { -// image1Size[0] = nX1; -// image1Size[1] = nY1; -// image2Size[0] = nX2; -// image2Size[1] = nY2; - customized_2DSIZE = true; if(m_DRTGeometryMetaInfo == NULL) { // todo @@ -334,12 +267,12 @@ void itkImageProcessor::SetCustom_2Dsize(int nX1, int nY1,int nX2,int nY2) ImageType3D::SizeType Size; Size [0] = nX1; Size [1] = nY1; - Size [1] = 1.; + Size [2] = 1.; m_DRTGeometryMetaInfo->SetDRT1Size(Size); Size [0] = nX2; Size [1] = nY2; - Size [1] = 1.; + Size [2] = 1.; m_DRTGeometryMetaInfo->SetDRT2Size(Size); //TODO UPDATE TO FOLLOW @@ -349,8 +282,6 @@ void itkImageProcessor::SetCustom_2Dsize(int nX1, int nY1,int nX2,int nY2) int itkImageProcessor::load3DSerie(const char * pcDirName) { - // std::cout<<"itkImageProcessor::load3DSerie"<SetImageIO(dicomIO); imageSeriesReader3D->SetFileNames(fileNames); @@ -440,7 +371,7 @@ int itkImageProcessor::load3DSerie(const char * pcDirName) if(!strcmp(sTmpString,ImageOrientationStrings[eImageOrientationType::HFS])){ m_PatOrientation = eImageOrientationType::HFS; } else if(!strcmp(sTmpString,ImageOrientationStrings[eImageOrientationType::FFS])){ - m_PatOrientation = eImageOrientationType::FFS; + m_PatOrientation = eImageOrientationType::FFS; } else { m_PatOrientation = eImageOrientationType::NotDefined; } @@ -460,6 +391,8 @@ int itkImageProcessor::load3DSerie(const char * pcDirName) m_VolumeSourceDupli->SetInputImage(caster3D->GetOutput()); m_VolumeSourceDupli->Update(); + caster3D = NULL; + } } @@ -507,10 +440,10 @@ int itkImageProcessor::load3DSerie(const char * pcDirName) /* initialise DRT meta */ m_DRTImage1MetaInfo = DRTImageMetaInformation::New(); m_DRTImage1MetaInfo->SetProjectionAngleLPS( - CalcProjectionAngleLPS( - m_CTMetaInfo->GetPatientOrientation(), - m_DRTGeometryMetaInfo->GetProjectionAngle1IEC()) - ); + CalcProjectionAngleLPS( + m_CTMetaInfo->GetPatientOrientation(), + m_DRTGeometryMetaInfo->GetProjectionAngle1IEC()) + ); m_DRTImage1MetaInfo->SetSCD( m_DRTGeometryMetaInfo->GetSCD()); @@ -519,17 +452,32 @@ int itkImageProcessor::load3DSerie(const char * pcDirName) m_DRTImage2MetaInfo = DRTImageMetaInformation::New(); m_DRTImage2MetaInfo->SetProjectionAngleLPS( - CalcProjectionAngleLPS( - m_CTMetaInfo->GetPatientOrientation(), - m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) - ); + CalcProjectionAngleLPS( + m_CTMetaInfo->GetPatientOrientation(), + m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) + ); m_DRTImage2MetaInfo->SetSCD( m_DRTGeometryMetaInfo->GetSCD()); + std::cout<< m_DRTGeometryMetaInfo->GetDRT1Size() <GetDRT1Spacing() <UpdateProjectionGeometryMeta(); + std::cout<< m_DRTGeometryMetaInfo->GetDRT1Size() <GetDRT1Spacing() <GetSize() <GetSpacing() <GetOrigin() <GetSize() <GetSpacing() <GetOrigin() <GetOutput(); - - - - /* TODO: Here calculate COV - * m_3DProjectionOriginLPSZero - * m_3DProjectionOriginLPS */ - - - return EXIT_SUCCESS; } itkImageProcessor::ImageType3D::PointType itkImageProcessor::CalculateDRTOrigin( - ImageType3D::SizeType m_ImageSize, - ImageType3D::SpacingType m_ImageResolution, - double dSCD + ImageType3D::SizeType m_ImageSize, + ImageType3D::SpacingType m_ImageResolution, + double dSCD ){ - ImageType3D::PointType Origin; - double o2Dx, o2Dy; + ImageType3D::PointType Origin; + double o2Dx, o2Dy; - o2Dx = ((double)m_ImageSize [0] - 1.) / 2.; - o2Dy = ((double)m_ImageSize [1] - 1.) / 2.; - // Compute the origin (in mm) of the 2D Image - Origin[0] = -m_ImageResolution [0] * o2Dx; - Origin[1] = -m_ImageResolution [1] * o2Dy; - Origin[2] = -dSCD ; + o2Dx = ((double)m_ImageSize [0] - 1.) / 2.; + o2Dy = ((double)m_ImageSize [1] - 1.) / 2.; + // Compute the origin (in mm) of the 2D Image + Origin[0] = -m_ImageResolution [0] * o2Dx; + Origin[1] = -m_ImageResolution [1] * o2Dy; + Origin[2] = -dSCD ; - return - Origin; + return + Origin; } itkImageProcessor::ImageType3D::PointType @@ -614,135 +553,49 @@ itkImageProcessor::CalculateProjectionCenterLPS( if(bZero) { return - m_ProjectionOriginLPS - Origin; + m_ProjectionOriginLPS - Origin; } else { return - m_ProjectionOriginLPS; + m_ProjectionOriginLPS; } } -//void itkImageProcessor::ApplyVolumeImportTransform(){ - -// //ImageType3D::PointType ImportOffset; - -// if(customized_ImportTransform && -// customized_RTCouchSetup && -// customized_RTIsocenter ){ - -// ImportOffset= -// this->CalcImportVolumeOffset(rtCouchOffset, -// m_CTMetaInfo->GetLPS2IECDirections(), -// rtIsocenterLPS, -// IEC2DCMMapT); -// std::cout<<"ImportOffset> "<GetOutput()->GetOrigin(); -// Origin = Origin - ImportOffset ; - -// ChangeInformationFilterType::Pointer -// m_3DInputChangeInformation = -// ChangeInformationFilterType::New(); -// m_3DInputChangeInformation->SetInput( -// m_VolumeSourceDupli->GetOutput()); -// m_3DInputChangeInformation->SetOutputOrigin(Origin ) ;//NewOrigin); -// m_3DInputChangeInformation->ChangeOriginOn(); -// m_3DInputChangeInformation->UpdateOutputInformation(); -// m_3DInputChangeInformation->Update(); - -// /* End of origin change */ - -// InternalImageType::Pointer m_Image = -// m_3DInputChangeInformation->GetOutput(); - -// // TODO: move this view origin to new meta variable. -// m_CTMetaInfo->SetOriginLPS(m_Image->GetOrigin()); - -// /* Calculate center of projection in patient coordinates */ - -// InternalImageType::DirectionType IECtoLPS_Directions; -// IECtoLPS_Directions = -// m_CTMetaInfo->GetLPS2IECDirections().GetTranspose(); -// // LPStoIEC_Directions.GetTranspose(); - -// if(customized_ProjectionCenter) -// { -// m_3DProjectionOriginLPS = -// IECtoLPS_Directions * m_DRTGeometryMetaInfo->GetProjectionCenter(); -// /* in longitudinal direction (Sup-Inf), we put it in the -// * middle of the volume */ -// m_3DProjectionOriginLPS[2] = m_CTMetaInfo->GetCOV()[2]; - -// } else { -// m_3DProjectionOriginLPS = m_CTMetaInfo->GetCOV(); -// } - - -// /* calculate also the same center with zero origin -// * This is required to bring back the DRT into LPS */ -// m_3DProjectionOriginLPSZero = -// m_3DProjectionOriginLPS - m_CTMetaInfo->GetOriginLPS() ; - -// /* We should know everything we need to calculate the -// * origin and direction of projection images */ - -// //std::cout<< "m_3DProjectionOriginLPS "<GetCOV() <GetProjectionOriginLPS()<GetProjectionOriginLPSZero()<GetProjectionOriginLPS()<GetProjectionOriginLPSZero()< itkImageProcessor::GetRTImportOffset() { + std::vector vOffset; + vOffset.clear(); + + vOffset.push_back( + m_CTMetaInfo->GetImportOffset()[0] + ); + + vOffset.push_back( + m_CTMetaInfo->GetImportOffset()[1] + ); + + vOffset.push_back( + m_CTMetaInfo->GetImportOffset()[2] + ); return - ImportOffset.GetDataPointer(); + vOffset; } const std::vector itkImageProcessor::GetLPStoProjectionGeoLPSOffset() { - std::vector vOffset; - + vOffset.clear(); vOffset.push_back( m_DRTImage1MetaInfo->GetLPStoProjectionGeoLPSOffset()[0] - ); + ); vOffset.push_back( m_DRTImage1MetaInfo->GetLPStoProjectionGeoLPSOffset()[1] - ); + ); vOffset.push_back( m_DRTImage1MetaInfo->GetLPStoProjectionGeoLPSOffset()[2] - ); - + ); return vOffset; - } double @@ -1052,57 +905,57 @@ double itkImageProcessor::GetLocalizerDisplayWindowLevel(int iImg) switch (iImg) { - case 0: + case 0: m_TImageMeta = m_TImage1MetaInfo; - break; + break; - case 1: + case 1: m_TImageMeta = m_TImage2MetaInfo; - break; + break; - default: - return 0; - break; + default: + return 0; + break; - } + } - if(m_TImageMeta == NULL){ - return 0; - } else { - return - m_TImageMeta->GetWLLevel(); - } + if(m_TImageMeta == NULL){ + return 0; + } else { + return + m_TImageMeta->GetWLLevel(); + } } - double itkImageProcessor::GetLocalizerDisplayWindowWidth(int iImg) +double itkImageProcessor::GetLocalizerDisplayWindowWidth(int iImg) { TopogramImageMetaInformation::Pointer m_TImageMeta; - switch (iImg) { + switch (iImg) { - case 0: - m_TImageMeta = m_TImage1MetaInfo; - break; + case 0: + m_TImageMeta = m_TImage1MetaInfo; + break; - case 1: - m_TImageMeta = m_TImage2MetaInfo; + case 1: + m_TImageMeta = m_TImage2MetaInfo; - break; + break; - default: - return 0; - break; + default: + return 0; + break; - } + } - if(m_TImageMeta == NULL){ - return 0; - } else { - return - m_TImageMeta->GetWLWindow(); - } + if(m_TImageMeta == NULL){ + return 0; + } else { + return + m_TImageMeta->GetWLWindow(); + } } @@ -1122,13 +975,13 @@ void itkImageProcessor::InitializeProjector() const double dtr = (atan(1.0) * 4.0) / 180.0; transform->SetRotation( dtr * RZero[0], - dtr * RZero[1], - dtr * RZero[2]); + dtr * RZero[1], + dtr * RZero[2]); if(m_DRTImage1MetaInfo->GetProjectionOriginLPSZero() != - m_DRTImage2MetaInfo->GetProjectionOriginLPSZero() ) { - // TODO - } + m_DRTImage2MetaInfo->GetProjectionOriginLPSZero() ) { + // TODO + } transform->SetCenter( m_DRTImage1MetaInfo->GetProjectionOriginLPSZero() ); @@ -1142,9 +995,7 @@ void itkImageProcessor::InitializeProjector() interpolator1->SetThreshold( m_DRTGeometryMetaInfo->GetIntensityThreshold() ); - interpolator1->SetTransform(transform); - interpolator1->Initialize(); // 2D Image 2 @@ -1157,7 +1008,6 @@ void itkImageProcessor::InitializeProjector() m_DRTGeometryMetaInfo->GetIntensityThreshold() ); interpolator2->SetTransform(transform); - interpolator2->Initialize(); } @@ -1229,9 +1079,9 @@ void itkImageProcessor::loadRTPlanInfo( void itkImageProcessor::UpdateProjectionGeometryMeta(){ - if(m_CTMetaInfo == NULL){ - //TODO. - } + if(m_CTMetaInfo == NULL){ + //TODO. + } m_DRTImage1MetaInfo->SetProjectionOriginLPS( CalculateProjectionCenterLPS( @@ -1339,59 +1189,15 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ To be very careful editing ... */ -itkImageProcessor::ImageType3D::PointType -itkImageProcessor::CalcDRTImageOrigin(itkImageProcessor::InternalImageType::Pointer m_Image, - itkImageProcessor::ImageType3D::PointType m_VolCOOV, - double dAngle, - InternalImageType::DirectionType stdDRT2LPS - ){ - - 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; - - - 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]; - - 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::InternalImageType::DirectionType DRT2LPS - = this->CalcDRTImageDirections(dAngle, stdDRT2LPS); - - ImageType3D::PointType NewOrigin = - m_VolCOOV + DRT2LPS * (m_Image ->GetOrigin() - DRTCOV); - - return - NewOrigin; -} itkImageProcessor::ImageType3D::PointType itkImageProcessor::CalcDRTImageOrigin( - ImageType3D::PointType m_DRTOrigin, - ImageType3D::PointType m_DRTCOV, - ImageType3D::PointType m_ProjectionCenterLPS, - double dAngle, - InternalImageType::DirectionType stdDRT2LPS - ){ + ImageType3D::PointType m_DRTOrigin, + ImageType3D::PointType m_DRTCOV, + ImageType3D::PointType m_ProjectionCenterLPS, + double dAngle, + InternalImageType::DirectionType stdDRT2LPS + ){ itkImageProcessor::InternalImageType::DirectionType DRT2LPS = this->CalcDRTImageDirections(dAngle, stdDRT2LPS); @@ -1423,8 +1229,8 @@ void itkImageProcessor::GetProjectionImages(){ // The transform is determined by the parameters and the rotation center. if(m_DRTImage1MetaInfo->GetProjectionOriginLPSZero() != - m_DRTImage2MetaInfo->GetProjectionOriginLPSZero() ) { - // TODO + m_DRTImage2MetaInfo->GetProjectionOriginLPSZero() ) { + // TODO } finalTransform->SetCenter( @@ -1448,6 +1254,7 @@ void itkImageProcessor::GetProjectionImages(){ interpolator1->Initialize(); resampleFilter1->SetInterpolator(interpolator1); + resampleFilter1->SetSize(m_DRTImage1MetaInfo->GetSize()); resampleFilter1->SetOutputSpacing(m_DRTImage1MetaInfo->GetSpacing()); resampleFilter1->SetOutputOrigin(m_DRTImage1MetaInfo->GetOrigin()); @@ -1471,11 +1278,9 @@ void itkImageProcessor::GetProjectionImages(){ resampleFilter2->SetOutputSpacing(m_DRTImage2MetaInfo->GetSpacing()); resampleFilter2->SetOutputOrigin(m_DRTImage2MetaInfo->GetOrigin()); - resampleFilter1->Update(); resampleFilter2->Update(); - /* here probably we should move everything back to dicom */ /* Projection images are in a fancy reference system... @@ -1511,6 +1316,7 @@ void itkImageProcessor::GetProjectionImages(){ filter2->ChangeOriginOn(); filter2->UpdateOutputInformation(); + filter1->Update(); filter2->Update(); @@ -1552,60 +1358,60 @@ itkImageProcessor::CalcImportVolumeOffset( void itkImageProcessor::Write2DImages(){ - using RescaleFilterType = itk::RescaleIntensityImageFilter; - using WriterType = itk::ImageFileWriter; + // using RescaleFilterType = itk::RescaleIntensityImageFilter; + // using WriterType = itk::ImageFileWriter; - if(image2D1Loaded) - { + // if(image2D1Loaded) + // { - // Rescale the intensity of the projection images to 0-255 for output. - RescaleFilterType::Pointer rescaler1 = RescaleFilterType::New(); - rescaler1->SetOutputMinimum(0); - rescaler1->SetOutputMaximum(255); - rescaler1->SetInput(m_PASourceDupli->GetOutput()); + // // Rescale the intensity of the projection images to 0-255 for output. + // RescaleFilterType::Pointer rescaler1 = RescaleFilterType::New(); + // rescaler1->SetOutputMinimum(0); + // rescaler1->SetOutputMaximum(255); + // rescaler1->SetInput(m_PASourceDupli->GetOutput()); - WriterType::Pointer writer1 = WriterType::New(); + // WriterType::Pointer writer1 = WriterType::New(); - writer1->SetFileName("/mnt/gdrive/Scratch/gfattori/2D1.tif"); - writer1->SetInput(rescaler1->GetOutput()); + // writer1->SetFileName("/mnt/gdrive/Scratch/gfattori/2D1.tif"); + // writer1->SetInput(rescaler1->GetOutput()); - try - { - std::cout << "Writing image 2D1" << std::endl; - writer1->Update(); - } - catch (itk::ExceptionObject & err) - { - std::cerr << "ERROR: ExceptionObject caught !" << std::endl; - std::cerr << err << std::endl; - } + // try + // { + // std::cout << "Writing image 2D1" << std::endl; + // writer1->Update(); + // } + // catch (itk::ExceptionObject & err) + // { + // std::cerr << "ERROR: ExceptionObject caught !" << std::endl; + // std::cerr << err << std::endl; + // } - } + // } - if(image2D2Loaded){ - // Rescale the intensity of the projection images to 0-255 for output. - RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New(); - rescaler2->SetOutputMinimum(0); - rescaler2->SetOutputMaximum(255); - rescaler2->SetInput(m_LATSourceDupli->GetOutput()); + // if(image2D2Loaded){ + // // Rescale the intensity of the projection images to 0-255 for output. + // RescaleFilterType::Pointer rescaler2 = RescaleFilterType::New(); + // rescaler2->SetOutputMinimum(0); + // rescaler2->SetOutputMaximum(255); + // rescaler2->SetInput(m_LATSourceDupli->GetOutput()); - WriterType::Pointer writer2 = WriterType::New(); + // WriterType::Pointer writer2 = WriterType::New(); - writer2->SetFileName("/mnt/gdrive/Scratch/gfattori/2D2.tif"); - writer2->SetInput(rescaler2->GetOutput()); + // writer2->SetFileName("/mnt/gdrive/Scratch/gfattori/2D2.tif"); + // writer2->SetInput(rescaler2->GetOutput()); - try - { - std::cout << "Writing image 2D2" << std::endl; - writer2->Update(); - } - catch (itk::ExceptionObject & err) - { - std::cerr << "ERROR: ExceptionObject caught !" << std::endl; - std::cerr << err << std::endl; - } - } + // try + // { + // std::cout << "Writing image 2D2" << std::endl; + // writer2->Update(); + // } + // catch (itk::ExceptionObject & err) + // { + // std::cerr << "ERROR: ExceptionObject caught !" << std::endl; + // std::cerr << err << std::endl; + // } + // } } @@ -1788,36 +1594,31 @@ vtkImageData* itkImageProcessor::GetProjection1VTK() vtkMatrix4x4 * itkImageProcessor::GetProjection1VTKTransform() { -// std::cout<< "Composed " <GetComposedTransform()->Print(std::cout); -// std::cout<< "Inverse " <GetInverseTransform()->Print(std::cout); + // std::cout<< "Composed " <GetComposedTransform()->Print(std::cout); + // std::cout<< "Inverse " <GetInverseTransform()->Print(std::cout); m_Projection1VTKTransform->Identity(); for(int iIdx = 0; iIdx<3 ; iIdx++){ for(int jIdx = 0; jIdx<3 ; jIdx++){ m_Projection1VTKTransform->SetElement(iIdx,jIdx, - interpolator1->GetInverseTransform()->GetMatrix().GetVnlMatrix()[iIdx][jIdx]); + interpolator1->GetInverseTransform()->GetMatrix().GetVnlMatrix()[iIdx][jIdx]); } } m_Projection1VTKTransform->SetElement(0,3, - interpolator1->GetInverseTransform()->GetTranslation()[0]); + interpolator1->GetInverseTransform()->GetTranslation()[0]); m_Projection1VTKTransform->SetElement(1,3, - interpolator1->GetInverseTransform()->GetTranslation()[1]); + interpolator1->GetInverseTransform()->GetTranslation()[1]); m_Projection1VTKTransform->SetElement(2,3, - interpolator1->GetInverseTransform()->GetTranslation()[2]); - - - //m_Projection1VTKTransform->Print(std::cout); + interpolator1->GetInverseTransform()->GetTranslation()[2]); return m_Projection1VTKTransform; - - } vtkMatrix4x4 * itkImageProcessor::GetProjection2VTKTransform() @@ -1828,18 +1629,18 @@ vtkMatrix4x4 * itkImageProcessor::GetProjection2VTKTransform() for(int iIdx = 0; iIdx<3 ; iIdx++){ for(int jIdx = 0; jIdx<3 ; jIdx++){ m_Projection2VTKTransform->SetElement(iIdx,jIdx, - interpolator2->GetInverseTransform()->GetMatrix().GetVnlMatrix()[iIdx][jIdx]); + interpolator2->GetInverseTransform()->GetMatrix().GetVnlMatrix()[iIdx][jIdx]); } } m_Projection2VTKTransform->SetElement(0,3, - interpolator2->GetInverseTransform()->GetTranslation()[0]); + interpolator2->GetInverseTransform()->GetTranslation()[0]); m_Projection2VTKTransform->SetElement(1,3, - interpolator2->GetInverseTransform()->GetTranslation()[1]); + interpolator2->GetInverseTransform()->GetTranslation()[1]); m_Projection2VTKTransform->SetElement(2,3, - interpolator2->GetInverseTransform()->GetTranslation()[2]); + interpolator2->GetInverseTransform()->GetTranslation()[2]); @@ -1977,7 +1778,6 @@ void itkImageProcessor::SetInitialTranslations(double dX, double dY, double dZ) IECtoLPS_Directions = m_CTMetaInfo->GetLPS2IECDirections().GetTranspose(); -// LPStoIEC_Directions.GetTranspose(); ImageType3D::PointType LPSTranslations = IECtoLPS_Directions * IECTranslations; @@ -1986,9 +1786,6 @@ void itkImageProcessor::SetInitialTranslations(double dX, double dY, double dZ) TZero[1]= LPSTranslations[1]; TZero[2]= LPSTranslations[2]; - - // std::cout<< "itkImageProcessor::SetInitialTranslations IEC : "< GetLPStoProjectionGeoLPSOffset(); /** Get Projection origin in LPS coordinates*/ - const double* GetRTImportOffset(); + const std::vector GetRTImportOffset(); /** Apply transform to CT volume. @@ -211,11 +209,6 @@ private: ImageType3D::PointType m_ProjectionCenter, bool bZero); - - ImageType3D::PointType ImportOffset; - ImageType3D :: PointType - ConversionOffset; - TransformType::Pointer transform; @@ -223,7 +216,6 @@ private: interpolator1, interpolator2; - ITKtoVTKFilterType::Pointer toVTK2D1, toVTK2D2, @@ -243,10 +235,6 @@ private: ChangeInformationFilterType::Pointer m_3DInputChangeInformationToZero; -// ImageType3D :: PointType -// m_3DProjectionOriginLPS, -// m_3DProjectionOriginLPSZero; - void UpdateProjectionGeometryMeta(); @@ -257,13 +245,6 @@ private: CalcDRTImageDirections(double angle, InternalImageType::DirectionType DRT2LPS); - ImageType3D::PointType - CalcDRTImageOrigin(InternalImageType::Pointer m_Image, - ImageType3D::PointType m_VolCOOV, - double dAngle, - InternalImageType::DirectionType stdDRT2LPS - ); - ImageType3D::PointType CalcDRTImageOrigin( ImageType3D::PointType m_DRTOrigin, @@ -271,7 +252,7 @@ private: ImageType3D::PointType m_ProjectionCenterLPS, double dAngle, InternalImageType::DirectionType stdDRT2LPS - ); + ); double @@ -286,49 +267,12 @@ private: ImageType3D::PointType rtIsocenterLPS, ImageType3D::PointType IEC2DCMMapT); -// double image1res[2], image2res[2]; -// int image1Size[2], image2Size[2]; - double TZero[3], RZero[3]; - bool - customized_2DRES, - customized_2DSIZE, - customized_RTIsocenter, - customized_RTCouchSetup, - customized_ImportTransform, - customized_ProjectionCenter; - - bool - image2D1Loaded, - image2D2Loaded, - image3DLoaded; - - - - /* Transform between IEC Support and - * IEC Scanner frame of reference */ -// ImageType3D::PointType -// IEC2DCMMapT, -// IEC2DCMMapR; - - /* RT Plan isocenter in LPS coordinates */ -// ImageType3D::PointType -// rtIsocenterLPS; - - /* RT Plan couch SETUP offset corresponding to - * the RT Plan isocenter */ -// ImageType3D::PointType -// rtCouchOffset; - - - InternalImageType::DirectionType Std_DRT2LPS; - - vtkMatrix4x4 * m_Projection1VTKTransform, * m_Projection2VTKTransform; @@ -355,6 +299,9 @@ private: RTGeometryMetaInformation::Pointer m_RTMetaInfo; + TransformMetaInformation::Pointer + m_TransformMetaInfo; + };