Commit 7f44f22c authored by Emmanuel Promayon's avatar Emmanuel Promayon

FIXED translation management

Uses cube face intersection to determine where (in %) is the
current translation
parent e33447ae
......@@ -47,29 +47,34 @@
namespace camitk {
// Useful debug macros for displaying homogeneous matrix and points
#define displayPoint(...) CAMITK_INFO_ALT(#__VA_ARGS__ + QString("\n[%1,%2,%3,%4]") \
.arg(__VA_ARGS__[0], 8, 'G', 5, ' ') \
.arg(__VA_ARGS__[1], 8, 'G', 5, ' ') \
.arg(__VA_ARGS__[2], 8, 'G', 5, ' ') \
.arg(__VA_ARGS__[3], 8, 'G', 5, ' '))
#define displayPoint(...) CAMITK_INFO_ALT(#__VA_ARGS__ + QString(" = [%1,%2,%3,%4]") \
.arg(__VA_ARGS__[0], 8, 'f', 4, ' ') \
.arg(__VA_ARGS__[1], 8, 'f', 4, ' ') \
.arg(__VA_ARGS__[2], 8, 'f', 4, ' ') \
.arg(__VA_ARGS__[3], 8, 'f', 4, ' '))
#define displayQVector3D(...) CAMITK_INFO_ALT(#__VA_ARGS__ + QString(" = (%1,%2,%3)") \
.arg(__VA_ARGS__.x(), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__.y(), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__.z(), 8, 'f', 4, ' '))
#define displayMatrix4x4(...) CAMITK_INFO_ALT(#__VA_ARGS__ + QString("\n[%1,%2,%3,%4]\n[%5,%6,%7,%8]\n[%9,%10,%11,%12]\n[%13,%14,%15,%16]") \
.arg(__VA_ARGS__->GetElement(0, 0), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(0, 1), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(0, 2), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(0, 3), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(1, 0), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(1, 1), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(1, 2), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(1, 3), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(2, 0), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(2, 1), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(2, 2), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(2, 3), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(3, 0), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(3, 1), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(3, 2), 8, 'G', 5, ' ') \
.arg(__VA_ARGS__->GetElement(3, 3), 8, 'G', 5, ' '))
.arg(__VA_ARGS__->GetElement(0, 0), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(0, 1), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(0, 2), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(0, 3), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(1, 0), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(1, 1), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(1, 2), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(1, 3), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(2, 0), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(2, 1), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(2, 2), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(2, 3), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(3, 0), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(3, 1), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(3, 2), 8, 'f', 4, ' ') \
.arg(__VA_ARGS__->GetElement(3, 3), 8, 'f', 4, ' '))
template<typename T>
vtkSmartPointer<vtkMatrix4x4> Multiply4x4(T a, T b) {
......@@ -82,18 +87,54 @@ template<typename T, typename... Args>
vtkSmartPointer<vtkMatrix4x4> Multiply4x4(T a, T b, Args... args) {
return Multiply4x4(a, Multiply4x4(b, args...));
}
QVector3D linePlaneIntersectionPoint(QVector3D lineVector, QVector3D linePoint, QVector3D planeNormal, QVector3D planePoint) {
double lineDPPlane = QVector3D::dotProduct(lineVector, planeNormal);
if (lineDPPlane < 1e-10) {
// line and plane are parallel (hopefully should never happen!)
return planePoint;
// bool linePlaneIntersection(Vector& contact, Vector ray, Vector rayOrigin, Vector normal, Vector coord) {
// // get d value
// float d = Dot(normal, coord);
// if (Dot(normal, ray) == 0) {
// return false; // No intersection, the line is parallel to the plane
// }
// // Compute the X value for the directed line ray intersecting the plane
// float x = (d - Dot(normal, rayOrigin)) / Dot(normal, ray);
// // output contact point
// *contact = rayOrigin + normalize(ray)*x; //Make sure your ray vector is normalized
// return true;
// }
// Vector3D diff = rayPoint - planePoint;
// double prod1 = diff.dot(planeNormal);
// double prod2 = rayVector.dot(planeNormal);
// double prod3 = prod1 / prod2;
// return rayPoint - rayVector * prod3;
bool linePlaneIntersectionPoint(QVector3D lineVector, QVector3D linePoint, QVector3D planeNormal, QVector3D planePoint, QVector3D& intersection) {
lineVector.normalize();
planeNormal.normalize();
// Let P(x,y,z) be the intersection point
// As the plane equation is:
// (P - planePoint) . planeNormal = 0 (. denotes dot product)
// and the line equation:
// P = linePoint + k * lineVector
// The equation linking both above is:
// (linePoint + k * lineVector - planePoint) . planeNormal = 0
// => k = - [ (linePoint - planePoint).planeNormal ] / (lineVector . planeNormal)
// if (lineVector . planeNormal) == 0.0 line is parallel to plane, this method should return false
float lDotN = QVector3D::dotProduct(lineVector, planeNormal);
displayQVector3D(lineVector)
displayQVector3D(linePoint)
displayQVector3D(planeNormal)
displayQVector3D(planePoint)
if (fabs(lDotN) < 1e-10) {
CAMITK_INFO_ALT(QString("line and plane are parallel, lDotN = %1").arg(lDotN));
// line and plane are parallel
return false;
}
else {
QVector3D planeToPoint = linePoint - planePoint;
double planeToPointDPPlane = QVector3D::dotProduct(planeToPoint, planeNormal);
double k = lineDPPlane / planeToPointDPPlane;
return linePoint - k * lineVector;
QVector3D u = linePoint - planePoint; // vector from plane point to the line point
float uDotN = QVector3D::dotProduct(u, planeNormal);
float k = - uDotN / lDotN;
intersection = linePoint + k * lineVector;
CAMITK_INFO_ALT(QString("lDotN = %1, uDotN = %2, k = %3").arg(lDotN).arg(uDotN).arg(k));
displayQVector3D(u)
displayQVector3D(intersection)
return true;
}
}
......@@ -131,9 +172,10 @@ void ArbitrarySingleImageComponent::resetTransform() {
getTransform()->Identity();
getTransform()->GetMatrix()->SetElement(0, 3, 0.0);
getTransform()->GetMatrix()->SetElement(1, 3, 0.0);
// TODO set to 0.5
getTransform()->GetMatrix()->SetElement(2, 3, dimensions[2] / 2.0 * spacing[2]);
getTransform()->Modified();
// TODO set to 0.5
// setTransformTranslation(0.0, 0.0, 0.5);
}
// -------------------- setTransformRotation --------------------
......@@ -147,14 +189,15 @@ void ArbitrarySingleImageComponent::setTransformRotation(double angleX, double a
vtkSmartPointer<vtkMatrix4x4> T_P2L = vtkSmartPointer<vtkMatrix4x4>::New();
T_P2L->DeepCopy(getTransform()->GetMatrix());
// this is only the rotation part of P2L, required for substracting current rotation in T_P2L
vtkSmartPointer<vtkMatrix4x4> T_P2L_rotation_only = vtkSmartPointer<vtkMatrix4x4>::New();
T_P2L_rotation_only->DeepCopy(T_P2L);
// this is only the rotation part of P2L, required for absolute rotation
T_P2L_rotation_only->SetElement(0, 3, 0.0);
T_P2L_rotation_only->SetElement(1, 3, 0.0);
T_P2L_rotation_only->SetElement(2, 3, 0.0);
T_P2L_rotation_only->SetElement(3, 3, 1.0);
// multiplying by T_P2L_rotation_only_inverse will remove the current rotation in T_P2L
vtkSmartPointer<vtkMatrix4x4> T_P2L_rotation_only_inverse = vtkSmartPointer<vtkMatrix4x4>::New();
vtkMatrix4x4::Invert(T_P2L_rotation_only, T_P2L_rotation_only_inverse);
......@@ -184,7 +227,7 @@ void ArbitrarySingleImageComponent::setTransformRotation(double angleX, double a
// Cumulative rotation = rotate(...)
// vtkSmartPointer<vtkMatrix4x4> checkRotation = Multiply4x4(T_P2C, R_P, T_C2P, T_P2L);
//checkRotation = by(by(by(by(T_P2C, R_P), T_P2L_rotation_only_inverse), T_C2P), T_P2L);
// Absolute rotation is as this:
vtkSmartPointer<vtkMatrix4x4> checkRotation = Multiply4x4(T_P2C, R_P, T_P2L_rotation_only_inverse, T_C2P, T_P2L);
if (checkCenter(checkRotation)) {
......@@ -200,69 +243,121 @@ void ArbitrarySingleImageComponent::setTransformRotation(double angleX, double a
// -------------------- setTransformTranslation --------------------
void ArbitrarySingleImageComponent::setTransformTranslation(double x, double y, double z) {
// Do nothing, my parent (the ImageComponent) should do that globally
// unless this is a arbitrary slice
if (sliceOrientation == Slice::ARBITRARY) {
// translation = set position to be at z% along the axis cMinus_P to cPlus_P
// ignore x and y
if (z < 0.0) {
z = 0.0;
}
else {
if (z > 1.0) {
z = 1.0;
}
}
CAMITK_INFO_ALT(QString("setTransformTranslation(%1)").arg(z))
updateTranslationExtremity();
QVector3D newCenter_P = cMinus_P + z * QVector3D(cPlus_P - cMinus_P);
displayQVector3D(newCenter_P)
if (!pointInsideVolume(newCenter_P)) {
CAMITK_INFO_ALT(QString("Can't translate to z=%1").arg(z))
return;
}
else {
vtkSmartPointer<vtkMatrix4x4> T_P2L = vtkSmartPointer<vtkMatrix4x4>::New();
T_P2L->DeepCopy(getTransform()->GetMatrix());
displayMatrix4x4(T_P2L)
vtkSmartPointer<vtkMatrix4x4> T_L2P = vtkSmartPointer<vtkMatrix4x4>::New();
vtkMatrix4x4::Invert(T_P2L, T_L2P);
// translation vector (x,y,z) is given in local coordinate system
double V_L[4] = {x, y, z, 1.0};
displayPoint(V_L)
double V_P[4];
T_P2L->MultiplyPoint(V_L, V_P);
displayPoint(V_P)
// translation of V expressed in P
vtkSmartPointer<vtkMatrix4x4> VTranslation_P = vtkSmartPointer<vtkMatrix4x4>::New();
VTranslation_P->Identity();
VTranslation_P->SetElement(0, 3, V_P[0]);
VTranslation_P->SetElement(1, 3, V_P[1]);
VTranslation_P->SetElement(2, 3, V_P[2]);
VTranslation_P->SetElement(3, 3, V_P[3]);
displayMatrix4x4(VTranslation_P)
vtkSmartPointer<vtkMatrix4x4> T_P2L_translation_only = vtkSmartPointer<vtkMatrix4x4>::New();
T_P2L_translation_only->Identity();
// this is only the translation part of L2P, required for absolute translation
T_P2L_translation_only->SetElement(0, 3, T_P2L->GetElement(0, 3));
T_P2L_translation_only->SetElement(1, 3, T_P2L->GetElement(1, 3));
T_P2L_translation_only->SetElement(2, 3, T_P2L->GetElement(2, 3));
T_P2L_translation_only->SetElement(3, 3, T_P2L->GetElement(3, 3));
displayMatrix4x4(T_P2L_translation_only)
vtkSmartPointer<vtkMatrix4x4> T_P2L_translation_only_inverse = vtkSmartPointer<vtkMatrix4x4>::New();
vtkMatrix4x4::Invert(T_P2L_translation_only, T_P2L_translation_only_inverse);
displayMatrix4x4(T_P2L_translation_only_inverse)
// cumulative translation
// vtkSmartPointer<vtkMatrix4x4> checkTranslation = vtkSmartPointer<vtkMatrix4x4>::New();
// vtkMatrix4x4::Multiply4x4(VTranslation_P,T_P2L,checkTranslation);
// displayMatrix4x4(checkTranslation)
vtkSmartPointer<vtkMatrix4x4> intermediateRotation1 = vtkSmartPointer<vtkMatrix4x4>::New();
vtkSmartPointer<vtkMatrix4x4> checkTranslation = vtkSmartPointer<vtkMatrix4x4>::New();
vtkMatrix4x4::Multiply4x4(VTranslation_P, T_P2L_translation_only_inverse, intermediateRotation1);
displayMatrix4x4(intermediateRotation1)
vtkMatrix4x4::Multiply4x4(intermediateRotation1, T_P2L, checkTranslation);
displayMatrix4x4(checkTranslation)
// move only if displacement will keep the center or any corners inside the box
if (checkCenter(checkTranslation)) {
CAMITK_INFO_ALT(QString("Ok to move to (%1,%2,%3)").arg(x).arg(y).arg(z))
getTransform()->GetMatrix()->DeepCopy(checkTranslation);
getTransform()->Modified();
}
else {
CAMITK_INFO_ALT(QString("Can't move to (%1,%2,%3)").arg(x).arg(y).arg(z))
}
double C_L[4] = {dimensions[0]* spacing[0] / 2.0 - spacing[0] / 2.0,
dimensions[1]* spacing[1] / 2.0 - spacing[1] / 2.0,
0.0,
1.0
};
double C_P[4];
T_P2L->MultiplyPoint(C_L, C_P);
QVector3D centerToOrigin_P(T_P2L->GetElement(0, 3) - C_P[0],
T_P2L->GetElement(1, 3) - C_P[1],
T_P2L->GetElement(2, 3) - C_P[2]);
displayQVector3D(centerToOrigin_P)
QVector3D newOrigin_P = newCenter_P + centerToOrigin_P;
displayQVector3D(newOrigin_P)
// Change the translation part only
T_P2L->SetElement(0, 3, newOrigin_P.x());
T_P2L->SetElement(1, 3, newOrigin_P.y());
T_P2L->SetElement(2, 3, newOrigin_P.z());
CAMITK_INFO_ALT(QString("Move to (%1,%2,%3)").arg(x).arg(y).arg(z))
getTransform()->GetMatrix()->DeepCopy(T_P2L);
getTransform()->Modified();
}
/*
// Do nothing, my parent (the ImageComponent) should do that globally
// unless this is a arbitrary slice
if (sliceOrientation == Slice::ARBITRARY) {
vtkSmartPointer<vtkMatrix4x4> T_P2L = vtkSmartPointer<vtkMatrix4x4>::New();
T_P2L->DeepCopy(getTransform()->GetMatrix());
displayMatrix4x4(T_P2L)
vtkSmartPointer<vtkMatrix4x4> T_L2P = vtkSmartPointer<vtkMatrix4x4>::New();
vtkMatrix4x4::Invert(T_P2L, T_L2P);
// translation vector (x,y,z) is given in local coordinate system
double V_L[4] = {x, y, z, 1.0};
displayPoint(V_L)
double V_P[4];
T_P2L->MultiplyPoint(V_L, V_P);
displayPoint(V_P)
// translation of V expressed in P
vtkSmartPointer<vtkMatrix4x4> VTranslation_P = vtkSmartPointer<vtkMatrix4x4>::New();
VTranslation_P->Identity();
VTranslation_P->SetElement(0, 3, V_P[0]);
VTranslation_P->SetElement(1, 3, V_P[1]);
VTranslation_P->SetElement(2, 3, V_P[2]);
VTranslation_P->SetElement(3, 3, V_P[3]);
displayMatrix4x4(VTranslation_P)
vtkSmartPointer<vtkMatrix4x4> T_P2L_translation_only = vtkSmartPointer<vtkMatrix4x4>::New();
T_P2L_translation_only->Identity();
// this is only the translation part of L2P, required for absolute translation
T_P2L_translation_only->SetElement(0, 3, T_P2L->GetElement(0, 3));
T_P2L_translation_only->SetElement(1, 3, T_P2L->GetElement(1, 3));
T_P2L_translation_only->SetElement(2, 3, T_P2L->GetElement(2, 3));
T_P2L_translation_only->SetElement(3, 3, T_P2L->GetElement(3, 3));
displayMatrix4x4(T_P2L_translation_only)
vtkSmartPointer<vtkMatrix4x4> T_P2L_translation_only_inverse = vtkSmartPointer<vtkMatrix4x4>::New();
vtkMatrix4x4::Invert(T_P2L_translation_only, T_P2L_translation_only_inverse);
displayMatrix4x4(T_P2L_translation_only_inverse)
// cumulative translation
// vtkSmartPointer<vtkMatrix4x4> checkTranslation = vtkSmartPointer<vtkMatrix4x4>::New();
// vtkMatrix4x4::Multiply4x4(VTranslation_P,T_P2L,checkTranslation);
// displayMatrix4x4(checkTranslation)
vtkSmartPointer<vtkMatrix4x4> intermediateRotation1 = vtkSmartPointer<vtkMatrix4x4>::New();
vtkSmartPointer<vtkMatrix4x4> checkTranslation = vtkSmartPointer<vtkMatrix4x4>::New();
vtkMatrix4x4::Multiply4x4(VTranslation_P, T_P2L_translation_only_inverse, intermediateRotation1);
displayMatrix4x4(intermediateRotation1)
vtkMatrix4x4::Multiply4x4(intermediateRotation1, T_P2L, checkTranslation);
displayMatrix4x4(checkTranslation)
// move only if displacement will keep the center or any corners inside the box
if (checkCenter(checkTranslation)) {
CAMITK_INFO_ALT(QString("Ok to move to (%1,%2,%3)").arg(x).arg(y).arg(z))
getTransform()->GetMatrix()->DeepCopy(checkTranslation);
getTransform()->Modified();
}
else {
CAMITK_INFO_ALT(QString("Can't move to (%1,%2,%3)").arg(x).arg(y).arg(z))
}
}
*/
}
// -------------------- updateTranslationExtremity --------------------
......@@ -279,50 +374,77 @@ void ArbitrarySingleImageComponent::updateTranslationExtremity() {
double C_P[4];
T_P2L->MultiplyPoint(C_L, C_P);
displayPoint(C_P)
// z vector in local coordinate system
double V_L[4] = { 0.0, 0.0, 1.0, 1.0};
double V_L[4] = { C_L[0], C_L[1], C_L[3] + 1.0, 1.0};
// z vector in global coordinate system
double V_P[4];
T_P2L->MultiplyPoint(V_L, V_P);
displayMatrix4x4(T_P2L)
displayPoint(V_P)
// compute intersection of line C_P + k V_P with the image
QVector3D lineVector(V_P[0], V_P[1], V_P[2]);
// compute intersection of line C_P + k vec(C_P, V_P) with the image
QVector3D lineVector(V_P[0] - C_P[0], V_P[1] - C_P[1], V_P[2] - C_P[2]);
QVector3D startPoint(C_P[0], C_P[1], C_P[2]);
QVector3D topCorner(0.0, 0.0, 0.0);
QVector3D bottomCorner(dimensions[0]*spacing[0], dimensions[1]*spacing[1], dimensions[2]*spacing[2]);
// check intersection to front/back plane
QVector3D intersection = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, 0.0, 1.0), topCorner);
if (intersection != topCorner && pointInsideVolume(intersection)) {
cPlus = intersection;
cMinus = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, 0.0, -1.0), bottomCorner);
double xCenter = dimensions[0] * spacing[0] / 2.0 - spacing[0] / 2.0;
double yCenter = dimensions[1] * spacing[1] / 2.0 - spacing[1] / 2.0;
double zCenter = dimensions[2] * spacing[2] / 2.0 - spacing[2] / 2.0;
// check intersection to front plane
double zMin = spacing[2] / 2.0;
QVector3D frontCenter(xCenter, yCenter, zMin);
QVector3D intersection;
CAMITK_INFO_ALT("----------- front")
bool intersect = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, 0.0, 1.0), frontCenter, intersection);
if (intersect && pointInsideVolume(intersection)) {
cMinus_P = intersection;
// back plane
double zMax = dimensions[2] * spacing[2] - spacing[2] / 2.0;
CAMITK_INFO_ALT("----------- back")
intersect = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, 0.0, -1.0), QVector3D(xCenter, yCenter, zMax), cPlus_P);
CAMITK_INFO_ALT("front/back")
}
else {
// check intersection to top/bottom plane
intersection = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, 1.0, 0.0), topCorner);
if (intersection != topCorner && pointInsideVolume(intersection)) {
cPlus = intersection;
cMinus = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, -1.0, 0.0), bottomCorner);
// check intersection to top plane
double yMin = spacing[1] / 2.0;
QVector3D topCenter(xCenter, yMin, zCenter);
CAMITK_INFO_ALT("----------- top")
intersect = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, 1.0, 0.0), topCenter, intersection);
displayQVector3D(intersection)
if (intersect && pointInsideVolume(intersection)) {
cMinus_P = intersection;
// bottom plane
double yMax = dimensions[1] * spacing[1] - spacing[1] / 2.0;
CAMITK_INFO_ALT("----------- bottom")
intersect = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(0.0, -1.0, 0.0), QVector3D(xCenter, yMax, zCenter), cPlus_P);
CAMITK_INFO_ALT("top/bottom")
}
else {
// intersection is with left/right plane
cPlus = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(1.0, 0.0, 0.0), topCorner);
cMinus = linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(-1.0, 0.0, 0.0), bottomCorner);
double xMin = spacing[0] / 2.0;
double xMax = dimensions[0] * spacing[0] - spacing[0] / 2.0;
CAMITK_INFO_ALT("----------- left")
linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(1.0, 0.0, 0.0), QVector3D(xMin, yCenter, zCenter), cMinus_P);
CAMITK_INFO_ALT("----------- right")
linePlaneIntersectionPoint(lineVector, startPoint, QVector3D(-1.0, 0.0, 0.0), QVector3D(xMax, yCenter, zCenter), cPlus_P);
CAMITK_INFO_ALT("left/right")
}
}
// TODO reorder cPlus/cMinus if needed
// TODO reorder cPlus_P/cMinus_P if needed
// C+ = C + k z
// if k < 0 inverse C+ and C-
CAMITK_INFO_ALT(QString("C- = (%1,%2,%3) C+ = (%4,%5,%6)").arg(cMinus.x()).arg(cMinus.y()).arg(cMinus.z()).arg(cPlus.x()).arg(cPlus.y()).arg(cPlus.z()))
// if k < 0 inverse C+ and C-
displayQVector3D(cMinus_P)
displayQVector3D(cPlus_P)
}
// -------------------- getTranslationInVolume --------------------
double ArbitrarySingleImageComponent::getTranslationInVolume() {
double ArbitrarySingleImageComponent::getTranslationInVolume() {
updateTranslationExtremity();
// compute center in P reference frame
double C_L[4] = {dimensions[0]* spacing[0] / 2.0 - spacing[0] / 2.0,
dimensions[1]* spacing[1] / 2.0 - spacing[1] / 2.0,
......@@ -336,12 +458,12 @@ double ArbitrarySingleImageComponent::getTranslationInVolume() {
double C_P[4];
T_P2L->MultiplyPoint(C_L, C_P);
QVector3D C(C_P[0], C_P[1], C_P[2]);
// C- = 0.0
// C+ = 1.0
// C-C = k C-C+
// k = C-C / C-C+
double k = QVector3D((C - cMinus) / (cPlus - cMinus)).length();
// k = C-C / C-C+
double k = QVector3D(C - cMinus_P).length() / QVector3D(cPlus_P - cMinus_P).length();
CAMITK_INFO_ALT(QString("k=%1").arg(k))
return k;
}
......@@ -368,9 +490,12 @@ bool ArbitrarySingleImageComponent::checkCenter(vtkSmartPointer<vtkMatrix4x4> tr
// -------------------- pointInsideVolume --------------------
bool ArbitrarySingleImageComponent::pointInsideVolume(QVector3D p) {
return (p.x() >= 0 && p.x() <= dimensions[0] * spacing[0]
&& p.y() >= 0 && p.y() <= dimensions[1] * spacing[1]
&& p.z() >= 0 && p.z() <= dimensions[2] * spacing[2]);
return (p.x() >= 0.0 /*+ spacing[0] / 2.0 */
&& p.x() <= dimensions[0] * spacing[0] - spacing[0] / 2.0
&& p.y() >= 0.0 /*+ spacing[1] / 2.0*/
&& p.y() <= dimensions[1] * spacing[1] - spacing[1] / 2.0
&& p.z() >= 0.0 /*+ spacing[2] / 2.0*/
&& p.z() <= dimensions[2] * spacing[2] - spacing[2] / 2.0);
}
......
......@@ -137,15 +137,17 @@ private:
/// point of intersection of the line passing at the center of the image
/// in the z+ direction, i.e. (0,0,1) in the local frame.
/// This is the extreme possible point to translate C in the z+ direction
/// This is the extreme possible point to translate C in the z+ direction.
/// This position is expressed in the parent frame coordinate system
/// (required to compute getTranslationInVolume and set new translation)
QVector3D cPlus;
QVector3D cPlus_P;
/// point of intersection of the line passing at the center of the image
/// in the z- direction, i.e. (0,0,-1) in the local frame.
/// This is the extreme possible point to translate C in the z- direction
/// This position is expressed in the parent frame coordinate system
/// (required to compute getTranslationInVolume and set new translation)
QVector3D cMinus;
QVector3D cMinus_P;
};
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment