Hi Carsten,
after reading a bit more about this, I'm aware why scale orientation is necessary in special cases. Our manipulation works now, too. Thanks for the support.
While playing around with Matrix::getTransform() we added the missing snuggle support. If you think it is worth to be commited, the patch is attached.
(Source may be optimized but this is more or less the original code from Graphics Gems IV).
Regards,
Michael
Gesendet: Freitag, 15. August 2014 um 19:45 Uhr
Von: "Michael Raab" <michael-r...@gmx.de>
An: opensg-users@lists.sourceforge.net
Betreff: Re: [Opensg-users] Problem with matrix decomposition
Von: "Michael Raab" <michael-r...@gmx.de>
An: opensg-users@lists.sourceforge.net
Betreff: Re: [Opensg-users] Problem with matrix decomposition
Hi Carsten,
Am 15.08.2014 16:51 schrieb Carsten Neumann <carsten.p.neum...@gmail.com>:
>
> Hello Michael,
>
> On 2014-08-15 14:09, Michael Raab wrote:
> > If I try to decompose 2 Matrices which are nearly the same I get a
> > result scale vector that has switch axis. Here's an example:
> > Matrix 1:
> > 1.040 -1.040 0.000 0.000
> > 0.707 0.707 0.000 0.000
> > 0.000 0.000 1.000 0.000
> > 0.000 0.000 0.000 1.000
> > Result 1:
> > Rotation 0, 0, 45
> > Scale 1.47031, 1, 1
> > Scale Orientation 0, 0, -0.382683, 0.92388
> > Matrix 2:
> > 1.046 -1.046 0.000 0.000
> > 0.707 0.707 0.000 0.000
> > 0.000 0.000 1.000 0.000
> > 0.000 0.000 0.000 1.000
> > Result 2:
> > Rotation: 0, 0, 45
> > Scale: 1, 1.47986, 1
> > Scale Orientation: 0, 0, 0.382683, 0.92388
> > Has someone an explanation for this?
>
> assuming that recombining this back into matrices yields the starting
> matrices (i.e. it is not a bug in the decomposition) this would be
> caused by a degree of freedom in the algorithm.
Does it really produce the initial Matrix? Haven't tested this...
> I can think of the following options:
> - figure out where this degree of freedom is and make choose the scale
> orientation more consistently.
Sorry I don't understand this..
> - can you use a ComponentTransform instead of a plain Transform and only
> modify the individual components. That would make the decomposed form
> the canonical representation and only build the matrix from those.
That is not an option as the matrix is the result of a sequence of transformations.
> - can you make assumptions about the matrix, e.g. is it composed from
> translation, rotation, and (uniform) scale only?
Unfortunately not. The only assumption I can think of is that the scale should be applied without any additional scale orientation.
Usually I think of transformation matrices as a combination of the 3 components translation, rotation and scale. Actually I do not know where the scale orientation is necessary or useful.
Michael
>
> Cheers,
> Carsten
>
> ------------------------------------------------------------------------------
> _______________________________________________
> Opensg-users mailing list
> Opensg-users@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/opensg-users
------------------------------------------------------------------------------
_______________________________________________
Opensg-users mailing list
Opensg-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/opensg-users
Am 15.08.2014 16:51 schrieb Carsten Neumann <carsten.p.neum...@gmail.com>:
>
> Hello Michael,
>
> On 2014-08-15 14:09, Michael Raab wrote:
> > If I try to decompose 2 Matrices which are nearly the same I get a
> > result scale vector that has switch axis. Here's an example:
> > Matrix 1:
> > 1.040 -1.040 0.000 0.000
> > 0.707 0.707 0.000 0.000
> > 0.000 0.000 1.000 0.000
> > 0.000 0.000 0.000 1.000
> > Result 1:
> > Rotation 0, 0, 45
> > Scale 1.47031, 1, 1
> > Scale Orientation 0, 0, -0.382683, 0.92388
> > Matrix 2:
> > 1.046 -1.046 0.000 0.000
> > 0.707 0.707 0.000 0.000
> > 0.000 0.000 1.000 0.000
> > 0.000 0.000 0.000 1.000
> > Result 2:
> > Rotation: 0, 0, 45
> > Scale: 1, 1.47986, 1
> > Scale Orientation: 0, 0, 0.382683, 0.92388
> > Has someone an explanation for this?
>
> assuming that recombining this back into matrices yields the starting
> matrices (i.e. it is not a bug in the decomposition) this would be
> caused by a degree of freedom in the algorithm.
Does it really produce the initial Matrix? Haven't tested this...
> I can think of the following options:
> - figure out where this degree of freedom is and make choose the scale
> orientation more consistently.
Sorry I don't understand this..
> - can you use a ComponentTransform instead of a plain Transform and only
> modify the individual components. That would make the decomposed form
> the canonical representation and only build the matrix from those.
That is not an option as the matrix is the result of a sequence of transformations.
> - can you make assumptions about the matrix, e.g. is it composed from
> translation, rotation, and (uniform) scale only?
Unfortunately not. The only assumption I can think of is that the scale should be applied without any additional scale orientation.
Usually I think of transformation matrices as a combination of the 3 components translation, rotation and scale. Actually I do not know where the scale orientation is necessary or useful.
Michael
>
> Cheers,
> Carsten
>
> ------------------------------------------------------------------------------
> _______________________________________________
> Opensg-users mailing list
> Opensg-users@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/opensg-users
------------------------------------------------------------------------------
_______________________________________________
Opensg-users mailing list
Opensg-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/opensg-users
Source/Base/Base/OSGGLFuncProtos.h | 1 + Source/Base/Base/OSGMatrix.h | 4 +- Source/Base/Base/OSGMatrix.inl | 192 +++++++++++++++++++++++++++++++++++++ Source/Base/Base/OSGVector.h | 4 +- Source/Base/Base/OSGVector.inl | 18 ++++ 5 files changed, 217 insertions(+), 2 deletions(-)
diff --git a/Source/Base/Base/OSGGLFuncProtos.h b/Source/Base/Base/OSGGLFuncProtos.h
index 750ea73..bd45043 100644
--- a/Source/Base/Base/OSGGLFuncProtos.h
+++ b/Source/Base/Base/OSGGLFuncProtos.h
@@ -679,6 +679,7 @@ typedef void (OSG_APIENTRY *osgGlBlendColorProc )(GLclampf,
GLclampf,
GLclampf);
typedef void (OSG_APIENTRY *osgGlBlendEquationProc )(GLenum);
+typedef void (OSG_APIENTRY *osgGlBlendEquationExtProc )(GLenum);
typedef void (OSG_APIENTRY *osgGlBlendFuncSeparateProc)(GLenum,
GLenum,
GLenum,
diff --git a/Source/Base/Base/OSGMatrix.h b/Source/Base/Base/OSGMatrix.h
index 99c4c57..703aedc 100644
--- a/Source/Base/Base/OSGMatrix.h
+++ b/Source/Base/Base/OSGMatrix.h
@@ -475,7 +475,9 @@ class TransformationMatrix
QuaternionType &r,
QuaternionType &so,
VectorType3f &s ) const;
-
+ QuaternionType snuggle(QuaternionType q,
+ VectorType3f &k ) const;
+
#ifdef __sgi
#pragma set woff 1424
#endif
diff --git a/Source/Base/Base/OSGMatrix.inl b/Source/Base/Base/OSGMatrix.inl
index 7f17365..1c31720 100644
--- a/Source/Base/Base/OSGMatrix.inl
+++ b/Source/Base/Base/OSGMatrix.inl
@@ -829,6 +829,9 @@ inline void
S.spectralDecompose(SO, s);
so.setValue(SO);
+
+ //QuaternionType p = snuggle(so, s);
+ //so.mult(p);
}
#ifdef __sgi
@@ -1383,6 +1386,195 @@ void TransformationMatrix<ValueTypeT>::getTransform(
scaleFactor *= flip;
}
+/******* Spectral Axis Adjustment *******/
+
+ /* Given a unit quaternion, q, and a scale vector, k, find a unit quaternion, p,
+ * which permutes the axes and turns freely in the plane of duplicate scale
+ * factors, such that q p has the largest possible w component, i.e. the
+ * smallest possible angle. Permutes k's components to go with q p instead of q.
+ * See Ken Shoemake and Tom Duff. Matrix Animation and Polar Decomposition.
+ * Proceedings of Graphics Interface 1992. Details on p. 262-263.
+ */
+template<class ValueTypeT> inline
+typename TransformationMatrix<ValueTypeT>::QuaternionType TransformationMatrix<ValueTypeT>::snuggle(
+ QuaternionType q,
+ VectorType3f& k) const
+{
+#define SQRTHALF (0.7071067811865475244)
+#define sgn(n,v) ((n)?-(v):(v))
+#define swap3(a,i,j) {a[3]=a[i]; a[i]=a[j]; a[j]=a[3];}
+#define cycle(a,p) if (p) {a[3]=a[0]; a[0]=a[1]; a[1]=a[2]; a[2]=a[3];} else {a[3]=a[2]; a[2]=a[1]; a[1]=a[0]; a[0]=a[3];}
+enum QuatPart {X, Y, Z, W};
+
+QuaternionType qxtoz(0,SQRTHALF,0,SQRTHALF);
+QuaternionType qytoz(SQRTHALF,0,0,SQRTHALF);
+QuaternionType qppmm(0.5, 0.5,-0.5,-0.5);
+QuaternionType qpppp(0.5, 0.5, 0.5, 0.5);
+QuaternionType qmpmm(-0.5, 0.5,-0.5,-0.5);
+QuaternionType qpppm(0.5, 0.5, 0.5,-0.5);
+QuaternionType q0001(0.0, 0.0, 0.0, 1.0);
+QuaternionType q1000(1.0, 0.0, 0.0, 0.0);
+
+ QuaternionType p = q0001;
+ double ka[4];
+ int i, turn = -1;
+ ka[X] = k.x(); ka[Y] = k.y(); ka[Z] = k.z();
+
+ if (ka[X]==ka[Y])
+ {
+ if (ka[X]==ka[Z])
+ turn = W;
+ else turn = Z;
+ }
+ else
+ {
+ if (ka[X]==ka[Z])
+ turn = Y;
+ else if (ka[Y]==ka[Z])
+ turn = X;
+ }
+ if (turn>=0)
+ {
+ QuaternionType qtoz, qp;
+ unsigned int win;
+ double mag[3], t;
+ switch (turn)
+ {
+ //default: return (Qt_Conj(q));
+ default: return q.conj();
+ //case X: q = Qt_Mul(q, qtoz = qxtoz); swap3(ka,X,Z) break;
+ case X: qtoz = qxtoz; q.mult(qtoz); swap3(ka,X,Z) break;
+ //case Y: q = Qt_Mul(q, qtoz = qytoz); swap3(ka,Y,Z) break;
+ case Y: qtoz = qytoz; q.mult(qtoz); swap3(ka,Y,Z) break;
+ case Z: qtoz = q0001; break;
+ }
+ //q = Qt_Conj(q);
+ q = q.conj();
+ mag[0] = (double)q.z()*q.z()+(double)q.w()*q.w()-0.5;
+ mag[1] = (double)q.x()*q.z()-(double)q.y()*q.w();
+ mag[2] = (double)q.y()*q.z()+(double)q.x()*q.w();
+
+ bool neg[3];
+ for (i=0; i<3; i++)
+ {
+ neg[i] = (mag[i]<0.0);
+ if (neg[i]) mag[i] = -mag[i];
+ }
+
+ if (mag[0]>mag[1])
+ {
+ if (mag[0]>mag[2])
+ win = 0;
+ else win = 2;
+ }
+ else
+ {
+ if (mag[1]>mag[2]) win = 1;
+ else win = 2;
+ }
+
+ switch (win)
+ {
+ case 0: if (neg[0]) p = q1000; else p = q0001; break;
+ case 1: if (neg[1]) p = qppmm; else p = qpppp; cycle(ka,0) break;
+ case 2: if (neg[2]) p = qmpmm; else p = qpppm; cycle(ka,1) break;
+ }
+
+ //qp = Qt_Mul(q, p);
+ qp = q;
+ qp.mult(p);
+ t = std::sqrt(mag[win]+0.5);
+ //p = Qt_Mul(p, Qt_(0.0,0.0,-qp.z/t,qp.w/t));
+ p.mult(QuaternionType(0.0,0.0,-qp.z()/t,qp.w()/t));
+ //p = Qt_Mul(qtoz, Qt_Conj(p));
+ QuaternionType qqq = qtoz;
+ qqq.mult(p.conj());
+ p = qqq;
+ }
+ else
+ {
+ double qa[4], pa[4];
+ unsigned int lo, hi;
+ bool par = false;
+ bool neg[4];
+ double all, big, two;
+ qa[0] = q.x(); qa[1] = q.y(); qa[2] = q.z(); qa[3] = q.w();
+ for (i=0; i<4; i++)
+ {
+ pa[i] = 0.0;
+ neg[i] = (qa[i]<0.0);
+ if (neg[i]) qa[i] = -qa[i];
+ par ^= neg[i];
+ }
+
+ // Find two largest components, indices in hi and lo
+ if (qa[0]>qa[1]) lo = 0;
+ else lo = 1;
+
+ if (qa[2]>qa[3]) hi = 2;
+ else hi = 3;
+
+ if (qa[lo]>qa[hi])
+ {
+ if (qa[lo^1]>qa[hi])
+ {
+ hi = lo; lo ^= 1;
+ }
+ else
+ {
+ hi ^= lo; lo ^= hi; hi ^= lo;
+ }
+ }
+ else
+ {
+ if (qa[hi^1]>qa[lo]) lo = hi^1;
+ }
+
+ all = (qa[0]+qa[1]+qa[2]+qa[3])*0.5;
+ two = (qa[hi]+qa[lo])*SQRTHALF;
+ big = qa[hi];
+ if (all>two)
+ {
+ if (all>big) //all
+ {
+ {int i; for (i=0; i<4; i++) pa[i] = sgn(neg[i], 0.5);}
+ cycle(ka,par);
+ }
+ else //big
+ {
+ pa[hi] = sgn(neg[hi],1.0);
+ }
+ }
+ else
+ {
+ if (two>big) //two
+ {
+ pa[hi] = sgn(neg[hi], SQRTHALF);
+ pa[lo] = sgn(neg[lo], SQRTHALF);
+ if (lo>hi)
+ {
+ hi ^= lo; lo ^= hi; hi ^= lo;
+ }
+ if (hi==W)
+ {
+ hi = "\001\002\000"[lo];
+ lo = 3-hi-lo;
+ }
+ swap3(ka,hi,lo);
+ }
+ else //big
+ {
+ pa[hi] = sgn(neg[hi],1.0);
+ }
+ }
+ p[0] = -pa[0]; p[1] = -pa[1]; p[2] = -pa[2]; p[3] = pa[3];
+ }
+
+ k[0] = ka[X]; k[1] = ka[Y]; k[2] = ka[Z];
+
+ return p;
+}
+
/*! \brief Factors a matrix m into 5 pieces: m = r s rt u t, where rt
means transpose of r, and r and u are rotations, s is a scale,
and t is a translation. Any projection information is returned
diff --git a/Source/Base/Base/OSGVector.h b/Source/Base/Base/OSGVector.h
index 86d747a..7ba09fd 100644
--- a/Source/Base/Base/OSGVector.h
+++ b/Source/Base/Base/OSGVector.h
@@ -505,7 +505,9 @@ class Point : public SelectVecStorage<ValueTypeT, SizeI>::Type
Point operator + (const VectorType &vec) const;
- Point operator - (const VectorType &vec) const;
+ Point operator - (const VectorType &vec) const;
+
+ Point operator + (const PointType &vec) const;
Point operator * (const ValueType val) const;
diff --git a/Source/Base/Base/OSGVector.inl b/Source/Base/Base/OSGVector.inl
index 02f4a78..4942851 100644
--- a/Source/Base/Base/OSGVector.inl
+++ b/Source/Base/Base/OSGVector.inl
@@ -1045,6 +1045,24 @@ Point<ValueTypeT, SizeI>
return returnValue;
}
+//! Component wise binary point addition operator
+
+template <class ValueTypeT,
+ UInt32 SizeI > inline
+ Point<ValueTypeT, SizeI>
+ Point<ValueTypeT, SizeI>::operator + (const PointType &pnt) const
+{
+ Point<ValueTypeT, SizeI> returnValue;
+
+ for(UInt32 i = 0; i < Self::_uiSize; i++)
+ {
+ returnValue[i] = Self::_values[i] + pnt[i];
+ }
+
+ return returnValue;
+}
+
+
//! Component wise binary scalar multiplication
template <class ValueTypeT,
------------------------------------------------------------------------------ Slashdot TV. Video for Nerds. Stuff that matters. http://tv.slashdot.org/
_______________________________________________ Opensg-users mailing list Opensg-users@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/opensg-users
