diff --git a/Code/Mantid/Framework/DataObjects/inc/MantidDataObjects/CoordTransformAffine.h b/Code/Mantid/Framework/DataObjects/inc/MantidDataObjects/CoordTransformAffine.h index 17e4d6f398597e84b9c1348df21edc4c82783028..05d9c3a09547a2ca0d82b52eb1204c3916fd33aa 100644 --- a/Code/Mantid/Framework/DataObjects/inc/MantidDataObjects/CoordTransformAffine.h +++ b/Code/Mantid/Framework/DataObjects/inc/MantidDataObjects/CoordTransformAffine.h @@ -43,6 +43,10 @@ public: const std::vector<Mantid::Kernel::VMD> &axes, const Mantid::Kernel::VMD &scaling); + void buildNonOrthogonal(const Mantid::Kernel::VMD &origin, + const std::vector<Mantid::Kernel::VMD> &axes, + const Mantid::Kernel::VMD &scaling); + virtual void apply(const coord_t *inputVector, coord_t *outVector) const; static CoordTransformAffine *combineTransformations(CoordTransform *first, diff --git a/Code/Mantid/Framework/DataObjects/src/CoordTransformAffine.cpp b/Code/Mantid/Framework/DataObjects/src/CoordTransformAffine.cpp index fcbe484f95ab17c881719ebfdc3c6069f93e753d..42d14ec96dc22fc6016767c0208d4cf43be78be0 100644 --- a/Code/Mantid/Framework/DataObjects/src/CoordTransformAffine.cpp +++ b/Code/Mantid/Framework/DataObjects/src/CoordTransformAffine.cpp @@ -206,6 +206,55 @@ void CoordTransformAffine::buildOrthogonal( copyRawMatrix(); } +void CoordTransformAffine::buildNonOrthogonal( + const Mantid::Kernel::VMD &origin, + const std::vector<Mantid::Kernel::VMD> &axes, + const Mantid::Kernel::VMD &scaling) { + if (origin.size() != inD) + throw std::runtime_error("CoordTransformAffine::buildNonOrthogonal(): the " + "origin must be in the dimensions of the input " + "workspace (length inD)."); + if (axes.size() != outD) + throw std::runtime_error("CoordTransformAffine::buildNonOrthogonal(): you " + "must give as many basis vectors as there are " + "dimensions in the output workspace."); + if (scaling.size() != outD) + throw std::runtime_error("CoordTransformAffine::buildNonOrthogonal(): the " + "size of the scaling vector must be the same as " + "the number of dimensions in the output " + "workspace."); + + // Start with identity + m_affineMatrix.identityMatrix(); + //A matrix is columns of basis vectors + Mantid::Kernel::Matrix<coord_t> A(inD,outD); + for(size_t i = 0; i < outD; i++){ + for(size_t j = 0; j < inD; j++){ + A[j][i]=axes[i][j]; + } + } + Mantid::Kernel::Matrix<coord_t> AT=A; + AT.Transpose(); + Mantid::Kernel::Matrix<coord_t> ATA=AT*A; + ATA.Invert(); + Mantid::Kernel::Matrix<coord_t> Ainv=ATA*AT; + Mantid::Kernel::Matrix<coord_t> offset(inD,1); + for(size_t j = 0; j < inD; j++){ + offset[j][0]=origin[j]; + } + Mantid::Kernel::Matrix<coord_t> outoffset=Ainv*offset; + + for(size_t i = 0; i < outD; i++){ + for(size_t j = 0; j < inD; j++){ + m_affineMatrix[i][j]=Ainv[i][j]*scaling[i]; + } + m_affineMatrix[i][inD]=-outoffset[i][0]*scaling[i]; + } + // Copy into the raw matrix (for speed) + copyRawMatrix(); +} + + //---------------------------------------------------------------------------------------------- /** Apply the coordinate transformation * diff --git a/Code/Mantid/Framework/Kernel/src/Matrix.cpp b/Code/Mantid/Framework/Kernel/src/Matrix.cpp index e1bb4e3e777f63f58fbf74c2d8bbbce3480d52d7..6a48da2f8c17c5812b0ea402ae04790ac165df13 100644 --- a/Code/Mantid/Framework/Kernel/src/Matrix.cpp +++ b/Code/Mantid/Framework/Kernel/src/Matrix.cpp @@ -965,6 +965,13 @@ T Matrix<T>::Invert() if (nx != ny && nx < 1) return 0; + if(nx==1) + { + T det=V[0][0]; + if(V[0][0]!=static_cast<T>(0.)) + V[0][0]=static_cast<T>(1.)/V[0][0]; + return det; + } int *indx = new int[nx]; // Set in lubcmp double *col = new double[nx]; diff --git a/Code/Mantid/Framework/Kernel/test/MatrixTest.h b/Code/Mantid/Framework/Kernel/test/MatrixTest.h index d27c2752545a5db01a97b8d5f1eaa373d09d14bc..af7d78e4973d9a10895dfea1ec6e309f9be6fcd4 100644 --- a/Code/Mantid/Framework/Kernel/test/MatrixTest.h +++ b/Code/Mantid/Framework/Kernel/test/MatrixTest.h @@ -55,6 +55,10 @@ public: A[2][2]=-7.0; TS_ASSERT_DELTA(A.Invert(),105.0,1e-5); + Matrix<double> B(1,1); + B[0][0]=2.; + TS_ASSERT_DELTA(B.Invert(),2,1e-5); + TS_ASSERT_DELTA(B[0][0],0.5,1e-5); } void testIdent() diff --git a/Code/Mantid/Framework/MDAlgorithms/src/SlicingAlgorithm.cpp b/Code/Mantid/Framework/MDAlgorithms/src/SlicingAlgorithm.cpp index 447eef142551061e9f428d8e2f12508d2e3afb09..82101d88f8572d9569306a05d3a1335a0becf061 100644 --- a/Code/Mantid/Framework/MDAlgorithms/src/SlicingAlgorithm.cpp +++ b/Code/Mantid/Framework/MDAlgorithms/src/SlicingAlgorithm.cpp @@ -394,15 +394,19 @@ void SlicingAlgorithm::createGeneralTransform() { DataObjects::CoordTransformAffine *ct = new DataObjects::CoordTransformAffine(inD, m_outD); // Note: the scaling makes the coordinate correspond to a bin index - ct->buildOrthogonal(m_inputMinPoint, this->m_bases, - VMD(this->m_binningScaling)); + //ct->buildOrthogonal(m_inputMinPoint, this->m_bases, + // VMD(this->m_binningScaling)); + ct->buildNonOrthogonal(m_inputMinPoint, this->m_bases, + VMD(this->m_binningScaling)/VMD(m_transformScaling)); this->m_transform = ct; // Transformation original->binned DataObjects::CoordTransformAffine *ctFrom = new DataObjects::CoordTransformAffine(inD, m_outD); - ctFrom->buildOrthogonal(m_translation, this->m_bases, - VMD(m_transformScaling)); + //ctFrom->buildOrthogonal(m_translation, this->m_bases, + // VMD(m_transformScaling)); + ctFrom->buildNonOrthogonal(m_translation, this->m_bases, + VMD(m_transformScaling)/VMD(m_transformScaling)); m_transformFromOriginal = ctFrom; // Validate