vtkContourTopogramProjectFilter works without transforms

This commit is contained in:
2022-02-15 15:20:41 +01:00
parent 867d2949a3
commit 3bc86087b6
5 changed files with 84 additions and 49 deletions

View File

@ -430,7 +430,7 @@ void itkImageProcessor::ApplyVolumeImportTransform(){
ImageType3D::PointType ImportOffset;
//ImageType3D::PointType ImportOffset;
if(customized_ImportTransform &&
customized_RTCouchSetup &&
@ -522,6 +522,9 @@ void itkImageProcessor::ApplyVolumeImportTransform(){
/* We should know everything we need to calculate the
* origin and direction of projection images */
//std::cout<< "m_3DProjectionOriginLPS "<<m_3DProjectionOriginLPS <<std::endl;
//std::cout<< "m_3DProjectionOriginLPSZero "<<m_3DProjectionOriginLPSZero <<std::endl;
if(true){
std::cout<<" -------- 3D IMAGE --------"<<std::endl;
@ -540,10 +543,23 @@ void itkImageProcessor::ApplyVolumeImportTransform(){
}
const double* itkImageProcessor::GetProjectionOriginLPS()
const double* itkImageProcessor::GetRTImportOffset()
{
return
m_3DProjectionOriginLPS.GetDataPointer();
ImportOffset.GetDataPointer();
}
const double* itkImageProcessor::GetLPStoProjectionGeoLPSOffset()
{
ConversionOffset[0] =
- m_3DProjectionOriginLPS[0] + m_3DProjectionOriginLPSZero[0];
ConversionOffset[1] =
- m_3DProjectionOriginLPS[1] + m_3DProjectionOriginLPSZero[1];
ConversionOffset[2] =
- m_3DProjectionOriginLPS[2] + m_3DProjectionOriginLPSZero[2];
return
ConversionOffset.GetDataPointer();
}
double
@ -1410,6 +1426,11 @@ vtkImageData* itkImageProcessor::GetProjection1VTK()
vtkMatrix4x4 * itkImageProcessor::GetProjection1VTKTransform()
{
// std::cout<< "Composed " <<std::endl;
// interpolator1->GetComposedTransform()->Print(std::cout);
// std::cout<< "Inverse " <<std::endl;
// interpolator1->GetInverseTransform()->Print(std::cout);
m_Projection1VTKTransform->Identity();
for(int iIdx = 0; iIdx<3 ; iIdx++){

View File

@ -111,7 +111,11 @@ public:
void InitializeProjector();
/** Get Projection origin in LPS coordinates*/
const double* GetProjectionOriginLPS();
const double* GetLPStoProjectionGeoLPSOffset();
/** Get Projection origin in LPS coordinates*/
const double* GetRTImportOffset();
/** Apply transform to CT volume.
* This moves the volume based on RT Plan info
@ -202,7 +206,9 @@ private:
using ITKtoVTKFilterType = itk::ImageToVTKImageFilter<OutputImageType>;
ImageType3D::PointType ImportOffset;
ImageType3D :: PointType
ConversionOffset;
TransformType::Pointer
transform;

View File

@ -168,6 +168,7 @@ public:
/** Get a pointer to the Inverse Transform. used to calculate the rays */
itkGetConstObjectMacro(InverseTransform, TransformType);
itkGetConstObjectMacro(ComposedTransform, TransformType);
/** Set and get the focal point to isocenter distance in mm */
itkSetMacro(FocalPointToIsocenterDistance, double);
@ -225,6 +226,7 @@ protected:
TransformPointer m_Transform; // Displacement of the volume
// Overall inverse transform used to calculate the ray position in the input space
TransformPointer m_InverseTransform;
TransformPointer m_ComposedTransform; // Composed transform
// The threshold above which voxels along the ray path are integrated
double m_Threshold;
@ -237,7 +239,6 @@ private:
TransformPointer m_GantryRotTransform; // Gantry rotation transform
TransformPointer m_CamShiftTransform; // Camera shift transform camRotTransform
TransformPointer m_CamRotTransform; // Camera rotation transform
TransformPointer m_ComposedTransform; // Composed transform
PointType m_SourcePoint; // Coordinate of the source in the standard Z projection geometry
PointType m_SourceWorld; // Coordinate of the source in the world coordinate system
};

View File

@ -34,6 +34,16 @@ vtkObjectFactoryNewMacro(vtkContourTopogramProjectionFilter);
vtkContourTopogramProjectionFilter::vtkContourTopogramProjectionFilter()
{
dProjectionLPSOffset[0]=0.;
dProjectionLPSOffset[1]=0.;
dProjectionLPSOffset[2]=0.;
dImportOffsetLPS[0] = 0.;
dImportOffsetLPS[1] = 0.;
dImportOffsetLPS[2] = 0.;
dSCD = -1.;
// by default process active point scalars
this->SetInputArrayToProcess(
0, 0, 0, vtkDataObject::FIELD_ASSOCIATION_POINTS, vtkDataSetAttributes::SCALARS);
@ -80,62 +90,44 @@ int vtkContourTopogramProjectionFilter::RequestData(
PrjPoints->SetNumberOfPoints( points->GetNumberOfPoints() );
double pp0[3], pp[3], pp2[3];
double pp0[3], pp1[3], pp2[3];
for(unsigned long int ii = 0; ii < points->GetNumberOfPoints() ; ii++){
//std::cout<< "-----" << std::endl;
pp0[0] = points->GetPoint(ii)[0] - dIso[0];
pp0[1] = points->GetPoint(ii)[1] - dIso[1];
pp0[2] = points->GetPoint(ii)[2] - dIso[2];
/** input points are in LPS, should be corrected for import offset and mapped onto
* the projection geometry origin */
pp0[0] = points->GetPoint(ii)[0] + dProjectionLPSOffset [0] - dImportOffsetLPS[0];
pp0[1] = points->GetPoint(ii)[1] + dProjectionLPSOffset [1] - dImportOffsetLPS[1];
pp0[2] = points->GetPoint(ii)[2] + dProjectionLPSOffset [2] - dImportOffsetLPS[2];
m_Transform->GetInverse()->TransformPoint(pp0,
pp);
/** LPS points in projection can now be mapped to the
* projection geometry using the inverse of the transform */
m_Transform->GetInverse()->TransformPoint(pp0, pp1);
/** coorrect x coordinate by the ratio of distance from panel */
if(pp[2] > 0) {
//pp[0] = pp[0] * 570/ pp[2];
} else if(pp[2] < 0){
//pp[0] = pp[0] / 570/ pp[2];
} else {
}
/** coorrect x coordinate by the ratio of distance from source
* Basic intercept theorem */
pp1[0] = pp1[0] * dSCD / abs(pp1[2]);
/** Go on the virtual detector */
pp[2] = -570.;
pp1[2] = -dSCD;
/** Points can be mapped back to LPS coordinates */
m_Transform->TransformPoint(pp1, pp2);
m_Transform->TransformPoint(pp,
pp2);
pp2[0] += dIso[0];
pp2[1] += dIso[1];
pp2[2] += dIso[2];
/** and back to 'real' LPS */
pp2[0] -= dProjectionLPSOffset [0];
pp2[1] -= dProjectionLPSOffset [1];
pp2[2] -= dProjectionLPSOffset [2];
PrjPoints->InsertPoint(ii,pp2);
// std::cout<< ": "<<ii <<std::endl;
// std::cout<< points->GetPoint(ii)[0] << " "
// << points->GetPoint(ii)[1] << " "
// << points->GetPoint(ii)[2]
// <<std::endl;
// std::cout<< PrjPoints ->GetPoint(ii)[0] << " "
// << PrjPoints ->GetPoint(ii)[1] << " "
// << PrjPoints ->GetPoint(ii)[2]
// <<std::endl;
// std::cout<< "-----" << std::endl;
}
// m_Transform->Print(std::cout);
//m_PrjPoly = vtkSmartPointer<vtkPolyData>::New();
output->SetPoints( PrjPoints );
//output->SetPoints(points);
output->BuildCells();
vtkIdType numPoints = points->GetNumberOfPoints();
if (psInput)
@ -164,9 +156,18 @@ int vtkContourTopogramProjectionFilter::RequestData(
}
void vtkContourTopogramProjectionFilter::SetProjectionCenter(const double * dP)
void vtkContourTopogramProjectionFilter::SetImportOffsetLPS(const double * dP)
{
memcpy(dIso,dP,3*sizeof(double));
memcpy(dImportOffsetLPS,dP,3*sizeof(double));
}
void vtkContourTopogramProjectionFilter::SetSCD(const double dVal){
dSCD = dVal;
}
void vtkContourTopogramProjectionFilter::SetLPStoProjectionLPSGeometryOffset(const double * dP)
{
memcpy(dProjectionLPSOffset,dP,3*sizeof(double));
}
void vtkContourTopogramProjectionFilter::SetTransformMatrix( vtkMatrix4x4 *userMatrix)

View File

@ -20,7 +20,9 @@ public:
void SetTransformMatrix( vtkMatrix4x4 *);
void SetProjectionCenter(const double * );
void SetLPStoProjectionLPSGeometryOffset (const double * );
void SetSCD(const double dVal);
void SetImportOffsetLPS(const double *);
protected:
vtkContourTopogramProjectionFilter();
@ -36,7 +38,11 @@ private:
void operator=(const vtkContourTopogramProjectionFilter&) = delete;
double dIso[3];
double
dProjectionLPSOffset[3],
dImportOffsetLPS[3],
dSCD;
vtkSmartPointer<vtkTransform> m_Transform;
vtkSmartPointer<vtkPolyData> m_PrjPoly;
};