Geo Cal review

2023-05-31 16:12:36 +02:00
parent 43a1d129dd
commit 74a294e06c
6 changed files with 49 additions and 56 deletions
--- a/reg23Topograms/itkDTRrecon/DRTMetaInformation.cpp
+++ b/reg23Topograms/itkDTRrecon/DRTMetaInformation.cpp
@ -94,7 +94,7 @@ DRTImageMetaInformation(){

    this->m_OriginLPS.Fill(0.);

-    this->m_PanelOffsetLPS.Fill(0.);
+    this->m_PanelOffsetPGeo = 0.;

 }

--- a/reg23Topograms/itkDTRrecon/DRTMetaInformation.h
+++ b/reg23Topograms/itkDTRrecon/DRTMetaInformation.h
@ -204,8 +204,8 @@ public:
    itkSetMacro(ProjectionOriginLPSZero,PointType);
    itkGetMacro(ProjectionOriginLPSZero,PointType);

-    itkSetMacro(PanelOffsetLPS,PointType);
-    itkGetMacro(PanelOffsetLPS,PointType);
+    itkSetMacro(PanelOffsetPGeo,double);
+    itkGetMacro(PanelOffsetPGeo,double);

    itkSetMacro(ProjectionAngleLPS,double);
    itkGetMacro(ProjectionAngleLPS,double);
@ -235,12 +235,12 @@ protected:
    //m_Origin,
    m_OriginLPS,
    m_ProjectionOriginLPS,
-    m_ProjectionOriginLPSZero,
-    m_PanelOffsetLPS;
+    m_ProjectionOriginLPSZero;

    double
    m_ProjectionAngleLPS,
-    m_SCD;
+    m_SCD,
+    m_PanelOffsetPGeo;

    DirectionType
        m_ImageDirectionsLPS;
--- a/reg23Topograms/itkDTRrecon/autoreg/itkMutualInformationTwoImageToOneImageMetric.hxx
+++ b/reg23Topograms/itkDTRrecon/autoreg/itkMutualInformationTwoImageToOneImageMetric.hxx
@ -148,7 +148,7 @@ MutualInformationTwoImageToOneImageMetric<TFixedImage, TMovingImage>::GetValue()
    //auto numberOfSamples = static_cast<unsigned int>(numberOfPixels * 0.50); //100%
    // since we have a ROI, then we should not set allPixels to TRUE.
    //metric1->UseAllPixelsOn();
-    std::cout << "m_NumberOfHistogramBins " << m_NumberOfHistogramBins << std::endl;
+    //std::cout << "m_NumberOfHistogramBins " << m_NumberOfHistogramBins << std::endl;
    metric1->SetNumberOfHistogramBins(m_NumberOfHistogramBins);

 //    InternalImageType::IndexType pIndex;
--- a/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp
+++ b/reg23Topograms/itkDTRrecon/itkImageProcessor.cpp
@ -674,6 +674,31 @@ const std::vector <double> itkImageProcessor::GetLPStoProjectionGeoLPSOffset()
            vOffset;
 }

+double
+itkImageProcessor::CalcPanelOffsetPGeo(
+        tPatOrientation pOrient,
+        double dPOffset){
+
+    double dVal = dPOffset;//0.;
+
+    switch (pOrient) {
+        case tPatOrientation :: HFS:
+
+        break;
+
+    case tPatOrientation :: FFS:
+        dVal *= -1;
+        break;
+
+    default:
+        break;
+
+    }
+
+    return
+            dVal;
+}
+
 double
 itkImageProcessor::CalcProjectionAngleLPS(
        tPatOrientation  pOrient,
@ -1266,7 +1291,7 @@ void itkImageProcessor::InitializeProjector()
            m_DRTImage2MetaInfo == NULL ||
            m_DRTGeometryMetaInfo == NULL ||
            m_TransformMetaInfo == NULL     ){
-        return;
+
        //TODO
        return;
    }
@ -1277,12 +1302,6 @@ void itkImageProcessor::InitializeProjector()
    ImageType3D::PointType ZeroPoint;
    ZeroPoint.Fill(0.);

-
-    //    InternalImageType::DirectionType IECtoLPS_Directions;
-
-    //    IECtoLPS_Directions =
-    //            m_CTMetaInfo->GetLPS2IECDirections().GetTranspose();
-
    TransformType::Pointer CurrTransform;

    CurrTransform=    CalculateInternalTransformV3(
@ -1316,8 +1335,7 @@ void itkImageProcessor::InitializeProjector()
    interpolator1->SetThreshold(
                m_DRTGeometryMetaInfo->GetIntensityThreshold() );
    interpolator1->SetPanelOffset(
-                m_DRTImage1MetaInfo->GetPanelOffsetLPS()[0]);
-//                m_DRTGeometryMetaInfo->GetPanel1Offset());
+                m_DRTImage1MetaInfo->GetPanelOffsetPGeo());

    interpolator1->SetTransform(transform1);
    interpolator1->Initialize();
@ -1354,8 +1372,7 @@ void itkImageProcessor::InitializeProjector()
    interpolator2->SetThreshold(
                m_DRTGeometryMetaInfo->GetIntensityThreshold() );
    interpolator2->SetPanelOffset(
-                m_DRTImage2MetaInfo->GetPanelOffsetLPS()[0]);
-//                m_DRTGeometryMetaInfo->GetPanel2Offset());
+                m_DRTImage2MetaInfo->GetPanelOffsetPGeo());

    interpolator2->SetTransform(transform2);
    interpolator2->Initialize();
@ -1630,24 +1647,17 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){
                m_DRTGeometryMetaInfo->GetDRT1Spacing()
                );

-
-    ImageType3D::PointType PanelOffsetIEC1;
-
-    PanelOffsetIEC1[0] = -m_DRTGeometryMetaInfo->GetPanel1Offset();
-    PanelOffsetIEC1[1] = 0.0;
-    PanelOffsetIEC1[2] = 0.0;
-
-    m_DRTImage1MetaInfo->SetPanelOffsetLPS(
-                m_CTMetaInfo->ConvertIECPointToLPSPoint(
-                                PanelOffsetIEC1));
-
+    m_DRTImage1MetaInfo->SetPanelOffsetPGeo(
+                this->CalcPanelOffsetPGeo(m_CTMetaInfo->GetPatientOrientation(),
+                                         -m_DRTGeometryMetaInfo->GetPanel1Offset())
+                );

    // This HAS TO be calculated from the nominal isocenter.
    m_DRTImage1MetaInfo->SetOriginLPS(
                CalcDRTImageOrigin(
                    m_DRTImage1MetaInfo->GetOrigin(),
                    m_DRTImage1MetaInfo->GetCOV(),
-                    NominalIsocenter_wZcorrectionLPS,//m_DRTImage1MetaInfo->GetProjectionOriginLPS(),
+                    NominalIsocenter_wZcorrectionLPS,
                    this->CalcProjectionAngleLPS(
                        m_CTMetaInfo->GetPatientOrientation(),
                        m_DRTGeometryMetaInfo->GetProjectionAngle1IEC() ),
@ -1735,23 +1745,15 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){
                m_DRTGeometryMetaInfo->GetDRT2Spacing()
                );

-
-    ImageType3D::PointType PanelOffsetIEC2;
-
-    PanelOffsetIEC2[0] = -m_DRTGeometryMetaInfo->GetPanel2Offset();
-    PanelOffsetIEC2[1] = 0.0;
-    PanelOffsetIEC2[2] = 0.0;
-
-    m_DRTImage2MetaInfo->SetPanelOffsetLPS(
-                m_CTMetaInfo->ConvertIECPointToLPSPoint(
-                                PanelOffsetIEC2));
-
+    m_DRTImage2MetaInfo->SetPanelOffsetPGeo(
+                this->CalcPanelOffsetPGeo(m_CTMetaInfo->GetPatientOrientation(),
+                                         -m_DRTGeometryMetaInfo->GetPanel2Offset()));

    m_DRTImage2MetaInfo->SetOriginLPS(
                CalcDRTImageOrigin(
                    m_DRTImage2MetaInfo->GetOrigin(),
                    m_DRTImage2MetaInfo->GetCOV(),
-                    NominalIsocenter_wZcorrectionLPS,//m_DRTImage2MetaInfo->GetProjectionOriginLPS(),
+                    NominalIsocenter_wZcorrectionLPS,
                    this->CalcProjectionAngleLPS(
                        m_CTMetaInfo->GetPatientOrientation(),
                        m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) ,
@ -1759,19 +1761,6 @@ void itkImageProcessor::UpdateProjectionGeometryMeta(){
                    ));


-    m_DRTImage2MetaInfo->SetOriginLPS(
-                CalcDRTImageOrigin(
-                    m_DRTImage2MetaInfo->GetOrigin(),
-                    m_DRTImage2MetaInfo->GetCOV(),
-                    NominalIsocenter_wZcorrectionLPS,//m_DRTImage2MetaInfo->GetProjectionOriginLPS(),
-                    this->CalcProjectionAngleLPS(
-                        m_CTMetaInfo->GetPatientOrientation(),
-                        m_DRTGeometryMetaInfo->GetProjectionAngle2IEC()) ,
-                    Std_DRT2LPS
-                    )
-
-                );
-
    m_DRTImage2MetaInfo->SetImageDirectionsLPS(
                CalcDRTImageDirections(
                    this->CalcProjectionAngleLPS(
--- a/reg23Topograms/itkDTRrecon/itkImageProcessor.h
+++ b/reg23Topograms/itkDTRrecon/itkImageProcessor.h
@ -352,6 +352,10 @@ private:
            tPatOrientation pOrient,
            double pAngleIEC);

+    double
+        CalcPanelOffsetPGeo(
+                tPatOrientation pOrient,
+                double dPOffset);
    /** Apply transform to CT volume.
     *  This moves the volume based on RT Plan info
     *  if those are available. */
--- a/reg23Topograms/itkDTRrecon/itkgSiddonJacobsRayCastInterpolateImageFunction.hxx
+++ b/reg23Topograms/itkDTRrecon/itkgSiddonJacobsRayCastInterpolateImageFunction.hxx
@ -158,8 +158,8 @@ gSiddonJacobsRayCastInterpolateImageFunction<TInputImage, TCoordRep>::Evaluate(c
  SlidingSourcePoint[0] = 0.;
  SlidingSourcePoint[1] = point[1];
  SlidingSourcePoint[2] = 0.;
-  PointType SourceWorld = m_InverseTransform->TransformPoint(SlidingSourcePoint);

+  PointType SourceWorld = m_InverseTransform->TransformPoint(SlidingSourcePoint);
  //std::cout<<"SourceWorld: "<<SourceWorld<<std::endl;

  // Get ths input pointers