From 3da8413f4f47d491c8affad4c4bae717070ed7a5 Mon Sep 17 00:00:00 2001 From: Fattori Giovanni Date: Fri, 24 Jun 2022 23:28:36 +0200 Subject: [PATCH] Load CT (3D) volume from database * If database contains a CT, that can be imported. * Full pipeline from DataView to itkImageProcessor. * No checks on CT selection in DataView. * No feedback of data loaded from RegistrationView to DataView. --- .../itkDTRrecon/itkImageProcessor.cpp | 308 +++++++++++------- .../itkDTRrecon/itkImageProcessor.h | 17 +- 2 files changed, 207 insertions(+), 118 deletions(-) diff --git a/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp b/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp index 3e678f2..49f55c2 100644 --- a/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp +++ b/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp @@ -284,12 +284,12 @@ void itkImageProcessor::SetCustom_ImportTransform(double dTx, } void itkImageProcessor::SetCustom_ProjectionCenterOffsetFixedAxes_IEC( - double dX1, - double dY1, - double dZ1, - double dX2, - double dY2, - double dZ2) + double dX1, + double dY1, + double dZ1, + double dX2, + double dY2, + double dZ2) { typedef itk::Point PointType; @@ -369,6 +369,82 @@ void itkImageProcessor::SetCustom_2Dsize(int nX1, int nY1,int nX2,int nY2) //TODO UPDATE TO FOLLOW } +int itkImageProcessor::load3DSerieFromFiles( std::vector v_fnames){ + + + tPatOrientation m_PatOrientation + = tPatOrientation::NotDefined; + + try + { + + using ImageIOType = itk::GDCMImageIO; + ImageIOType::Pointer dicomIO = ImageIOType::New(); + ImageSeriesReaderType3D::Pointer + imageSeriesReader3D = ImageSeriesReaderType3D::New(); + + imageSeriesReader3D->SetImageIO(dicomIO); + imageSeriesReader3D->SetFileNames(v_fnames); + imageSeriesReader3D->SetNumberOfWorkUnits(8); + imageSeriesReader3D->ForceOrthogonalDirectionOff(); // properly read CTs with gantry tilt + + imageSeriesReader3D->Update(); + + if(v_fnames.size()>0){ + + gdcm::Reader R; + R.SetFileName(v_fnames.at(0).c_str()); + if(!R.Read()) + { + cerr<<"ERROR: cannot read file: "<SetInput(imageSeriesReader3D->GetOutput()); + caster3D->Update(); + + m_VolumeSourceDupli->SetInputImage(caster3D->GetOutput()); + m_VolumeSourceDupli->Update(); + + caster3D = NULL; + + + } + catch (itk::ExceptionObject & ex) + { + std::cout << ex << std::endl; + return EXIT_FAILURE; + } + + InternalImageType::Pointer m_InputImage = + m_VolumeSourceDupli->GetOutput(); + return + this->fill3DVolumeMeta(m_InputImage, m_PatOrientation); + +} int itkImageProcessor::load3DSerie(const char * pcDirName) { @@ -493,9 +569,17 @@ int itkImageProcessor::load3DSerie(const char * pcDirName) return EXIT_FAILURE; } - InternalImageType::Pointer m_InputImage = m_VolumeSourceDupli->GetOutput(); + return + this->fill3DVolumeMeta(m_InputImage, m_PatOrientation); + +} + +/** Fill Meta after 3D volume load */ +int itkImageProcessor::fill3DVolumeMeta( + InternalImageType::Pointer m_InputImage, + tPatOrientation m_PatOrientation){ if(m_CTMetaInfo != NULL){ @@ -597,9 +681,9 @@ int itkImageProcessor::load3DSerie(const char * pcDirName) image3DIn = m_3DInputChangeInformationToZero->GetOutput(); - return EXIT_SUCCESS; + } const std::vector itkImageProcessor::GetRTImportOffset() @@ -1176,7 +1260,7 @@ void itkImageProcessor::InitializeProjector() //AAAAAA TODOOOOO CERCA ImageType3D::PointType pFakeIsoLPS = - m_DRTImage1MetaInfo->GetProjectionOriginLPS(); + m_DRTImage1MetaInfo->GetProjectionOriginLPS(); //ImageType3D::PointType pFakeIsoLPS = // m_CTMetaInfo->GetProjectionOriginLPS( @@ -1249,7 +1333,7 @@ void itkImageProcessor::InitializeProjector() if(m_RTMetaInfo == NULL) { ImageType3D::PointType pFakeIsoLPS = - m_DRTImage2MetaInfo->GetProjectionOriginLPS(); + m_DRTImage2MetaInfo->GetProjectionOriginLPS(); //ImageType3D::PointType pFakeIsoLPS = // m_CTMetaInfo->GetProjectionOriginLPS( @@ -1317,17 +1401,17 @@ void itkImageProcessor::InitializeProjector() -// // 2D Image 2 -// interpolator2->SetProjectionAngle( -// dtr * -// m_DRTImage2MetaInfo->GetProjectionAngleLPS() ); -// interpolator2->SetFocalPointToIsocenterDistance( -// m_DRTImage2MetaInfo->GetSCD()); -// interpolator2->SetThreshold( -// m_DRTGeometryMetaInfo->GetIntensityThreshold() -// ); -// interpolator2->SetTransform(transform); -// interpolator2->Initialize(); + // // 2D Image 2 + // interpolator2->SetProjectionAngle( + // dtr * + // m_DRTImage2MetaInfo->GetProjectionAngleLPS() ); + // interpolator2->SetFocalPointToIsocenterDistance( + // m_DRTImage2MetaInfo->GetSCD()); + // interpolator2->SetThreshold( + // m_DRTGeometryMetaInfo->GetIntensityThreshold() + // ); + // interpolator2->SetTransform(transform); + // interpolator2->Initialize(); @@ -1518,10 +1602,10 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ ImageType3D::PointType NominalIsocenterZero_wZcorrectionLPS; - NominalIsocenterZero_wZcorrectionLPS = - m_CTMetaInfo->GetProjectionOriginLPSZero( - m_DRTGeometryMetaInfo->GetProjectionCenter() - ); + NominalIsocenterZero_wZcorrectionLPS = + m_CTMetaInfo->GetProjectionOriginLPSZero( + m_DRTGeometryMetaInfo->GetProjectionCenter() + ); std::cout<< "///////////////// PGEOM META BEG ///////////////" <GetProjectionOriginLPS() ); - std::cout<<"Bounds1: " + std::cout<<"Bounds1: " <SetSizeWithBounds( vBounds.data(), @@ -1593,8 +1677,8 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ m_DRTImage1MetaInfo->GetCOV(), NominalIsocenter_wZcorrectionLPS,//m_DRTImage1MetaInfo->GetProjectionOriginLPS(), this->CalcProjectionAngleLPS( - m_CTMetaInfo->GetPatientOrientation(), - m_DRTGeometryMetaInfo->GetProjectionAngle1IEC() ), + m_CTMetaInfo->GetPatientOrientation(), + m_DRTGeometryMetaInfo->GetProjectionAngle1IEC() ), //m_DRTGeometryMetaInfo->GetProjectionAngle1IEC(),// m_DRTImage1MetaInfo->GetProjectionAngleLPS(), Std_DRT2LPS ) @@ -1605,9 +1689,9 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ m_DRTImage1MetaInfo->SetImageDirectionsLPS( CalcDRTImageDirections( this->CalcProjectionAngleLPS( - m_CTMetaInfo->GetPatientOrientation(), - m_DRTGeometryMetaInfo->GetProjectionAngle1IEC()), - //m_DRTGeometryMetaInfo->GetProjectionAngle1IEC(),//m_DRTImage1MetaInfo->GetProjectionAngleLPS(), + m_CTMetaInfo->GetPatientOrientation(), + m_DRTGeometryMetaInfo->GetProjectionAngle1IEC()), + //m_DRTGeometryMetaInfo->GetProjectionAngle1IEC(),//m_DRTImage1MetaInfo->GetProjectionAngleLPS(), Std_DRT2LPS) ); @@ -1648,15 +1732,15 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ m_DRTImage2MetaInfo->SetProjectionOriginLPS( CalibratedIsocenterLPS); -// m_CTMetaInfo->GetProjectionOriginLPS( -// m_DRTGeometryMetaInfo->GetProjectionCenter2()) -// ); + // m_CTMetaInfo->GetProjectionOriginLPS( + // m_DRTGeometryMetaInfo->GetProjectionCenter2()) + // ); m_DRTImage2MetaInfo->SetProjectionOriginLPSZero( CalibratedIsocenterZeroLPS); -// m_CTMetaInfo->GetProjectionOriginLPSZero( -// m_DRTGeometryMetaInfo->GetProjectionCenter2()) -// ); + // m_CTMetaInfo->GetProjectionOriginLPSZero( + // m_DRTGeometryMetaInfo->GetProjectionCenter2()) + // ); vBounds.clear(); vBounds = m_CTMetaInfo->CalculateRegions( @@ -1664,12 +1748,12 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ ); std::cout<<"Bounds2: " - <SetSizeWithBounds( vBounds.data(), @@ -1683,9 +1767,9 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ m_DRTImage2MetaInfo->GetCOV(), NominalIsocenter_wZcorrectionLPS,//m_DRTImage2MetaInfo->GetProjectionOriginLPS(), this->CalcProjectionAngleLPS( - m_CTMetaInfo->GetPatientOrientation(), - m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) , - //m_DRTGeometryMetaInfo->GetProjectionAngle2IEC(),//m_DRTImage2MetaInfo->GetProjectionAngleLPS(), + m_CTMetaInfo->GetPatientOrientation(), + m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) , + //m_DRTGeometryMetaInfo->GetProjectionAngle2IEC(),//m_DRTImage2MetaInfo->GetProjectionAngleLPS(), Std_DRT2LPS ) ); @@ -1693,8 +1777,8 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){ m_DRTImage2MetaInfo->SetImageDirectionsLPS( CalcDRTImageDirections( this->CalcProjectionAngleLPS( - m_CTMetaInfo->GetPatientOrientation(), - m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) , + m_CTMetaInfo->GetPatientOrientation(), + m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) , //m_DRTGeometryMetaInfo->GetProjectionAngle2IEC(),//m_DRTImage2MetaInfo->GetProjectionAngleLPS(), Std_DRT2LPS) ); @@ -1773,8 +1857,8 @@ void itkImageProcessor::GetProjectionImages(){ ImageType3D::PointType pFakeIsoLPS = m_DRTImage1MetaInfo->GetProjectionOriginLPS(); - //m_CTMetaInfo->GetProjectionOriginLPS( - // m_DRTGeometryMetaInfo->GetProjectionCenter()); + //m_CTMetaInfo->GetProjectionOriginLPS( + // m_DRTGeometryMetaInfo->GetProjectionCenter()); std::cout<<"pFakeIsoLPS: "<GetAngleZ() ); - std::cout<< "itkImageProcessor::GetProjectionImages" <SetCenter( m_DRTImage1MetaInfo->GetProjectionOriginLPSZero()); @@ -1851,65 +1935,65 @@ void itkImageProcessor::GetProjectionImages(){ //CurrTransform->SetComputeZYX(true); //CurrTransform->SetIdentity(); - if(m_RTMetaInfo == NULL) - { + if(m_RTMetaInfo == NULL) + { - ImageType3D::PointType pFakeIsoLPS = - m_DRTImage2MetaInfo->GetProjectionOriginLPS(); - //m_CTMetaInfo->GetProjectionOriginLPS( - // m_DRTGeometryMetaInfo->GetProjectionCenter()); + ImageType3D::PointType pFakeIsoLPS = + m_DRTImage2MetaInfo->GetProjectionOriginLPS(); + //m_CTMetaInfo->GetProjectionOriginLPS( + // m_DRTGeometryMetaInfo->GetProjectionCenter()); - std::cout<<"pFakeIsoLPS: "<GetTranslations(), - m_TransformMetaInfo->GetRotations(), - m_DRTImage2MetaInfo->GetProjectionOriginLPS(), - pFakeIsoLPS, - ZeroPoint, - IECtoLPS_Directions - ); - - std::cout<<"pFakeIsoLPS: "<GetIECS2IECScannerR(), - m_TransformMetaInfo->GetTranslations(), - m_TransformMetaInfo->GetRotations(), - m_DRTImage2MetaInfo->GetProjectionOriginLPS(), - m_RTMetaInfo->GetIsocenterLPS() - m_CTMetaInfo->GetImportOffset(), - m_CTMetaInfo->GetImportOffset(), - IECtoLPS_Directions - ); - - } - - - - transform2->SetComputeZYX(true); - transform2->SetIdentity(); - - transform2->SetTranslation( - CurrTransform->GetTranslation()); - transform2->SetRotation( - CurrTransform->GetAngleX(), - CurrTransform->GetAngleY(), - CurrTransform->GetAngleZ() + CurrTransform = + CalculateInternalTransform( + ZeroPoint, + ZeroPoint, + m_TransformMetaInfo->GetTranslations(), + m_TransformMetaInfo->GetRotations(), + m_DRTImage2MetaInfo->GetProjectionOriginLPS(), + pFakeIsoLPS, + ZeroPoint, + IECtoLPS_Directions ); + std::cout<<"pFakeIsoLPS: "<SetCenter( - m_DRTImage2MetaInfo->GetProjectionOriginLPSZero()); + } else { + + CurrTransform = + CalculateInternalTransform( + ZeroPoint , + m_DRTGeometryMetaInfo->GetIECS2IECScannerR(), + m_TransformMetaInfo->GetTranslations(), + m_TransformMetaInfo->GetRotations(), + m_DRTImage2MetaInfo->GetProjectionOriginLPS(), + m_RTMetaInfo->GetIsocenterLPS() - m_CTMetaInfo->GetImportOffset(), + m_CTMetaInfo->GetImportOffset(), + IECtoLPS_Directions + ); + + } - //transform2 ->Print(std::cout); + + transform2->SetComputeZYX(true); + transform2->SetIdentity(); + + transform2->SetTranslation( + CurrTransform->GetTranslation()); + transform2->SetRotation( + CurrTransform->GetAngleX(), + CurrTransform->GetAngleY(), + CurrTransform->GetAngleZ() + ); + + + transform2 ->SetCenter( + m_DRTImage2MetaInfo->GetProjectionOriginLPSZero()); + + + //transform2 ->Print(std::cout); @@ -1928,10 +2012,10 @@ void itkImageProcessor::GetProjectionImages(){ imageDRT2In = filter2->GetOutput(); std::cout<< "itkImageProcessor::GetProjectionImages END" <GetInverseTransform()->GetMatrix()<GetInverseTransform()->GetTranslation()<GetInverseTransform()->GetMatrix()<GetInverseTransform()->GetTranslation()< ); + int load2D(const char *); void loadRTPlanInfo(double, double, double, double, double ,double); @@ -151,13 +153,13 @@ public: protected: /** Various pixel types */ - using InternalPixelType = float; - using PixelType2D = short; - using PixelType3D = short; - using OutputPixelType = unsigned char; + using InternalPixelType = float; + using PixelType2D = short; + using PixelType3D = short; + using OutputPixelType = unsigned char; - using ImageType3D = itk::Image; - using InternalImageType = itk::Image; + using ImageType3D = itk::Image; + using InternalImageType = itk::Image; itkImageProcessor(); virtual ~itkImageProcessor(); @@ -167,6 +169,9 @@ private: itkImageProcessor(const Self&); //purposely not implemented void operator=(const Self&); //purposely not implemented + /** Fill Meta after 3D volume load */ + int fill3DVolumeMeta(InternalImageType::Pointer, + tPatOrientation); /** Image types */ using ImageType2D = itk::Image;