root/Math/CramerInversionSym.icc

0001 // @(#)root/smatrix:$Id$
0002 // Authors: L. Moneta    2005
0003
0004
0005 /**********************************************************************
0006  *                                                                    *
0007  * Copyright (c) 2005 , LCG ROOT MathLib Team                         *
0008  *                                                                    *
0009  *                                                                    *
0010  **********************************************************************/
0011 //
0012 // Cramer optimized inversion for symmetric matrices up to size 5x5.
0013 // Code from ROOT TMatrixDCramerInv which originates from CLHEP
0014 // (original author Mark Fischler)
0015 //
0016 // Modified by L. Moneta 22/03/07: specialize only until 5x5 (before was up to 6x6)
0017 //  tests show that on 64 machines (like on slc4) it is faster the general method
0018 //
0019
0020 #ifndef ROOT_Math_CramerInversionSym_icc
0021 #define ROOT_Math_CramerInversionSym_icc
0022
0023 #ifndef ROOT_Math_Dinv
0024 #error "Do not use CramerInversionSym.icc directly. #include \"Math/Dinv.h\" instead."
0025 #endif // ROOT_Math_Dinv
0026
0027 #include <cmath>
0028
0029
0030 namespace ROOT {
0031
0032   namespace Math {
0033
0034
0035
0036
0037 //==============================================================================
0038 /**
0039    Inversion for a 3x3 symmetric matrix
0040  */
0041 template <class T>
0042 bool FastInverter<3>::Dinv(MatRepSym<T,3> & rhs) {
0043
0044   typedef T  Scalar;
0045
0046   // check matrix sizes ??
0047
0048
0049   const Scalar c00 = rhs[4] * rhs[8] - rhs[5] * rhs[5];
0050   const Scalar c01 = rhs[5] * rhs[2] - rhs[1] * rhs[8];
0051   const Scalar c02 = rhs[1] * rhs[5] - rhs[4] * rhs[2];
0052   const Scalar c11 = rhs[8] * rhs[0] - rhs[2] * rhs[2];
0053   const Scalar c12 = rhs[2] * rhs[1] - rhs[5] * rhs[0];
0054   const Scalar c22 = rhs[0] * rhs[4] - rhs[1] * rhs[1];
0055
0056   const Scalar t0  = std::abs(rhs[0]);
0057   const Scalar t1  = std::abs(rhs[1]);
0058   const Scalar t2  = std::abs(rhs[2]);
0059
0060   Scalar det;
0061   Scalar tmp;
0062
0063   if (t0 >= t1) {
0064     if (t2 >= t0) {
0065       tmp = rhs[2];
0066       det = c12*c01-c11*c02;
0067     } else {
0068       tmp = rhs[0];
0069       det = c11*c22-c12*c12;
0070     }
0071   } else if (t2 >= t1) {
0072     tmp = rhs[2];
0073     det = c12*c01-c11*c02;
0074   } else {
0075     tmp = rhs[1];
0076     det = c02*c12-c01*c22;
0077   }
0078
0079   if ( det == 0 || tmp == 0)
0080     return false;
0081
0082   Scalar s = tmp/det;
0083 //   if (determ)
0084 //     *determ = 1./s;
0085
0086   rhs[0] = s*c00;
0087   rhs[1] = s*c01;
0088   rhs[2] = s*c02;
0089   rhs[4] = s*c11;
0090   rhs[5] = s*c12;
0091   rhs[8] = s*c22;
0092
0093   return true;
0094 }
0095
0096
0097 //==============================================================================
0098 // Inversion for 4x4 matrices
0099 //==============================================================================
0100
0101 // SFij are indices for a 4x4 symmetric matrix.
0102
0103 #define SF00 0
0104 #define SF01 1
0105 #define SF02 2
0106 #define SF03 3
0107
0108 #define SF10 1
0109 #define SF11 5
0110 #define SF12 6
0111 #define SF13 7
0112
0113 #define SF20 2
0114 #define SF21 6
0115 #define SF22 10
0116 #define SF23 11
0117
0118 #define SF30 3
0119 #define SF31 7
0120 #define SF32 11
0121 #define SF33 15
0122
0123
0124 /**
0125    Inversion for a 4x4 symmetric matrix
0126  */
0127 template <class T>
0128 bool FastInverter<4>::Dinv(MatRepSym<T,4> & rhs) {
0129
0130   typedef T  Scalar;
0131
0132
0133   // Find all NECESSARY 2x2 dets:  (14 of them)
0134
0135   const Scalar mDet2_12_01 = rhs[SF10]*rhs[SF21] - rhs[SF11]*rhs[SF20];
0136   const Scalar mDet2_12_02 = rhs[SF10]*rhs[SF22] - rhs[SF12]*rhs[SF20];
0137   const Scalar mDet2_12_12 = rhs[SF11]*rhs[SF22] - rhs[SF12]*rhs[SF21];
0138   const Scalar mDet2_13_01 = rhs[SF10]*rhs[SF31] - rhs[SF11]*rhs[SF30];
0139   const Scalar mDet2_13_02 = rhs[SF10]*rhs[SF32] - rhs[SF12]*rhs[SF30];
0140   const Scalar mDet2_13_03 = rhs[SF10]*rhs[SF33] - rhs[SF13]*rhs[SF30];
0141   const Scalar mDet2_13_12 = rhs[SF11]*rhs[SF32] - rhs[SF12]*rhs[SF31];
0142   const Scalar mDet2_13_13 = rhs[SF11]*rhs[SF33] - rhs[SF13]*rhs[SF31];
0143   const Scalar mDet2_23_01 = rhs[SF20]*rhs[SF31] - rhs[SF21]*rhs[SF30];
0144   const Scalar mDet2_23_02 = rhs[SF20]*rhs[SF32] - rhs[SF22]*rhs[SF30];
0145   const Scalar mDet2_23_03 = rhs[SF20]*rhs[SF33] - rhs[SF23]*rhs[SF30];
0146   const Scalar mDet2_23_12 = rhs[SF21]*rhs[SF32] - rhs[SF22]*rhs[SF31];
0147   const Scalar mDet2_23_13 = rhs[SF21]*rhs[SF33] - rhs[SF23]*rhs[SF31];
0148   const Scalar mDet2_23_23 = rhs[SF22]*rhs[SF33] - rhs[SF23]*rhs[SF32];
0149
0150   // SFind all NECESSSFRY 3x3 dets:   (10 of them)
0151
0152   const Scalar mDet3_012_012 = rhs[SF00]*mDet2_12_12 - rhs[SF01]*mDet2_12_02
0153                                + rhs[SF02]*mDet2_12_01;
0154   const Scalar mDet3_013_012 = rhs[SF00]*mDet2_13_12 - rhs[SF01]*mDet2_13_02
0155                                + rhs[SF02]*mDet2_13_01;
0156   const Scalar mDet3_013_013 = rhs[SF00]*mDet2_13_13 - rhs[SF01]*mDet2_13_03
0157                                + rhs[SF03]*mDet2_13_01;
0158   const Scalar mDet3_023_012 = rhs[SF00]*mDet2_23_12 - rhs[SF01]*mDet2_23_02
0159                                + rhs[SF02]*mDet2_23_01;
0160   const Scalar mDet3_023_013 = rhs[SF00]*mDet2_23_13 - rhs[SF01]*mDet2_23_03
0161                                + rhs[SF03]*mDet2_23_01;
0162   const Scalar mDet3_023_023 = rhs[SF00]*mDet2_23_23 - rhs[SF02]*mDet2_23_03
0163                                + rhs[SF03]*mDet2_23_02;
0164   const Scalar mDet3_123_012 = rhs[SF10]*mDet2_23_12 - rhs[SF11]*mDet2_23_02
0165                                + rhs[SF12]*mDet2_23_01;
0166   const Scalar mDet3_123_013 = rhs[SF10]*mDet2_23_13 - rhs[SF11]*mDet2_23_03
0167                                + rhs[SF13]*mDet2_23_01;
0168   const Scalar mDet3_123_023 = rhs[SF10]*mDet2_23_23 - rhs[SF12]*mDet2_23_03
0169                                + rhs[SF13]*mDet2_23_02;
0170   const Scalar mDet3_123_123 = rhs[SF11]*mDet2_23_23 - rhs[SF12]*mDet2_23_13
0171                                + rhs[SF13]*mDet2_23_12;
0172
0173   // Find the 4x4 det:
0174
0175   const Scalar det = rhs[SF00]*mDet3_123_123 - rhs[SF01]*mDet3_123_023
0176                        + rhs[SF02]*mDet3_123_013 - rhs[SF03]*mDet3_123_012;
0177
0178 //   if (determ)
0179 //     *determ = det;
0180
0181   if ( det == 0 )
0182     return false;
0183
0184   const Scalar oneOverDet = 1.0f / det;
0185   const Scalar mn1OverDet = - oneOverDet;
0186
0187   rhs[SF00] =  mDet3_123_123 * oneOverDet;
0188   rhs[SF01] =  mDet3_123_023 * mn1OverDet;
0189   rhs[SF02] =  mDet3_123_013 * oneOverDet;
0190   rhs[SF03] =  mDet3_123_012 * mn1OverDet;
0191
0192   rhs[SF11] =  mDet3_023_023 * oneOverDet;
0193   rhs[SF12] =  mDet3_023_013 * mn1OverDet;
0194   rhs[SF13] =  mDet3_023_012 * oneOverDet;
0195
0196   rhs[SF22] =  mDet3_013_013 * oneOverDet;
0197   rhs[SF23] =  mDet3_013_012 * mn1OverDet;
0198
0199   rhs[SF33] =  mDet3_012_012 * oneOverDet;
0200
0201   return true;
0202 }
0203
0204
0205 //==============================================================================
0206 // Inversion for 5x5 matrices
0207 //==============================================================================
0208
0209 // Mij are indices for a 5x5 matrix.
0210
0211 #define SM00 0
0212 #define SM01 1
0213 #define SM02 2
0214 #define SM03 3
0215 #define SM04 4
0216
0217 #define SM10 1
0218 #define SM11 6
0219 #define SM12 7
0220 #define SM13 8
0221 #define SM14 9
0222
0223 #define SM20 2
0224 #define SM21 7
0225 #define SM22 12
0226 #define SM23 13
0227 #define SM24 14
0228
0229 #define SM30 3
0230 #define SM31 8
0231 #define SM32 13
0232 #define SM33 18
0233 #define SM34 19
0234
0235 #define SM40 4
0236 #define SM41 9
0237 #define SM42 14
0238 #define SM43 19
0239 #define SM44 24
0240
0241 /**
0242    Inversion for a 5x5 symmetric matrix
0243  */
0244 template <class T>
0245 bool FastInverter<5>::Dinv(MatRepSym<T,5> & rhs) {
0246
0247   typedef T  Scalar;
0248
0249   // Find all NECESSARY 2x2 dets:  (25 of them)
0250
0251   const Scalar mDet2_23_01 = rhs[SM20]*rhs[SM31] - rhs[SM21]*rhs[SM30];
0252   const Scalar mDet2_23_02 = rhs[SM20]*rhs[SM32] - rhs[SM22]*rhs[SM30];
0253   const Scalar mDet2_23_03 = rhs[SM20]*rhs[SM33] - rhs[SM23]*rhs[SM30];
0254   const Scalar mDet2_23_12 = rhs[SM21]*rhs[SM32] - rhs[SM22]*rhs[SM31];
0255   const Scalar mDet2_23_13 = rhs[SM21]*rhs[SM33] - rhs[SM23]*rhs[SM31];
0256   const Scalar mDet2_23_23 = rhs[SM22]*rhs[SM33] - rhs[SM23]*rhs[SM32];
0257   const Scalar mDet2_24_01 = rhs[SM20]*rhs[SM41] - rhs[SM21]*rhs[SM40];
0258   const Scalar mDet2_24_02 = rhs[SM20]*rhs[SM42] - rhs[SM22]*rhs[SM40];
0259   const Scalar mDet2_24_03 = rhs[SM20]*rhs[SM43] - rhs[SM23]*rhs[SM40];
0260   const Scalar mDet2_24_04 = rhs[SM20]*rhs[SM44] - rhs[SM24]*rhs[SM40];
0261   const Scalar mDet2_24_12 = rhs[SM21]*rhs[SM42] - rhs[SM22]*rhs[SM41];
0262   const Scalar mDet2_24_13 = rhs[SM21]*rhs[SM43] - rhs[SM23]*rhs[SM41];
0263   const Scalar mDet2_24_14 = rhs[SM21]*rhs[SM44] - rhs[SM24]*rhs[SM41];
0264   const Scalar mDet2_24_23 = rhs[SM22]*rhs[SM43] - rhs[SM23]*rhs[SM42];
0265   const Scalar mDet2_24_24 = rhs[SM22]*rhs[SM44] - rhs[SM24]*rhs[SM42];
0266   const Scalar mDet2_34_01 = rhs[SM30]*rhs[SM41] - rhs[SM31]*rhs[SM40];
0267   const Scalar mDet2_34_02 = rhs[SM30]*rhs[SM42] - rhs[SM32]*rhs[SM40];
0268   const Scalar mDet2_34_03 = rhs[SM30]*rhs[SM43] - rhs[SM33]*rhs[SM40];
0269   const Scalar mDet2_34_04 = rhs[SM30]*rhs[SM44] - rhs[SM34]*rhs[SM40];
0270   const Scalar mDet2_34_12 = rhs[SM31]*rhs[SM42] - rhs[SM32]*rhs[SM41];
0271   const Scalar mDet2_34_13 = rhs[SM31]*rhs[SM43] - rhs[SM33]*rhs[SM41];
0272   const Scalar mDet2_34_14 = rhs[SM31]*rhs[SM44] - rhs[SM34]*rhs[SM41];
0273   const Scalar mDet2_34_23 = rhs[SM32]*rhs[SM43] - rhs[SM33]*rhs[SM42];
0274   const Scalar mDet2_34_24 = rhs[SM32]*rhs[SM44] - rhs[SM34]*rhs[SM42];
0275   const Scalar mDet2_34_34 = rhs[SM33]*rhs[SM44] - rhs[SM34]*rhs[SM43];
0276
0277   // Find all NECESSARY 3x3 dets:   (30 of them)
0278
0279   const Scalar mDet3_123_012 = rhs[SM10]*mDet2_23_12 - rhs[SM11]*mDet2_23_02 + rhs[SM12]*mDet2_23_01;
0280   const Scalar mDet3_123_013 = rhs[SM10]*mDet2_23_13 - rhs[SM11]*mDet2_23_03 + rhs[SM13]*mDet2_23_01;
0281   const Scalar mDet3_123_023 = rhs[SM10]*mDet2_23_23 - rhs[SM12]*mDet2_23_03 + rhs[SM13]*mDet2_23_02;
0282   const Scalar mDet3_123_123 = rhs[SM11]*mDet2_23_23 - rhs[SM12]*mDet2_23_13 + rhs[SM13]*mDet2_23_12;
0283   const Scalar mDet3_124_012 = rhs[SM10]*mDet2_24_12 - rhs[SM11]*mDet2_24_02 + rhs[SM12]*mDet2_24_01;
0284   const Scalar mDet3_124_013 = rhs[SM10]*mDet2_24_13 - rhs[SM11]*mDet2_24_03 + rhs[SM13]*mDet2_24_01;
0285   const Scalar mDet3_124_014 = rhs[SM10]*mDet2_24_14 - rhs[SM11]*mDet2_24_04 + rhs[SM14]*mDet2_24_01;
0286   const Scalar mDet3_124_023 = rhs[SM10]*mDet2_24_23 - rhs[SM12]*mDet2_24_03 + rhs[SM13]*mDet2_24_02;
0287   const Scalar mDet3_124_024 = rhs[SM10]*mDet2_24_24 - rhs[SM12]*mDet2_24_04 + rhs[SM14]*mDet2_24_02;
0288   const Scalar mDet3_124_123 = rhs[SM11]*mDet2_24_23 - rhs[SM12]*mDet2_24_13 + rhs[SM13]*mDet2_24_12;
0289   const Scalar mDet3_124_124 = rhs[SM11]*mDet2_24_24 - rhs[SM12]*mDet2_24_14 + rhs[SM14]*mDet2_24_12;
0290   const Scalar mDet3_134_012 = rhs[SM10]*mDet2_34_12 - rhs[SM11]*mDet2_34_02 + rhs[SM12]*mDet2_34_01;
0291   const Scalar mDet3_134_013 = rhs[SM10]*mDet2_34_13 - rhs[SM11]*mDet2_34_03 + rhs[SM13]*mDet2_34_01;
0292   const Scalar mDet3_134_014 = rhs[SM10]*mDet2_34_14 - rhs[SM11]*mDet2_34_04 + rhs[SM14]*mDet2_34_01;
0293   const Scalar mDet3_134_023 = rhs[SM10]*mDet2_34_23 - rhs[SM12]*mDet2_34_03 + rhs[SM13]*mDet2_34_02;
0294   const Scalar mDet3_134_024 = rhs[SM10]*mDet2_34_24 - rhs[SM12]*mDet2_34_04 + rhs[SM14]*mDet2_34_02;
0295   const Scalar mDet3_134_034 = rhs[SM10]*mDet2_34_34 - rhs[SM13]*mDet2_34_04 + rhs[SM14]*mDet2_34_03;
0296   const Scalar mDet3_134_123 = rhs[SM11]*mDet2_34_23 - rhs[SM12]*mDet2_34_13 + rhs[SM13]*mDet2_34_12;
0297   const Scalar mDet3_134_124 = rhs[SM11]*mDet2_34_24 - rhs[SM12]*mDet2_34_14 + rhs[SM14]*mDet2_34_12;
0298   const Scalar mDet3_134_134 = rhs[SM11]*mDet2_34_34 - rhs[SM13]*mDet2_34_14 + rhs[SM14]*mDet2_34_13;
0299   const Scalar mDet3_234_012 = rhs[SM20]*mDet2_34_12 - rhs[SM21]*mDet2_34_02 + rhs[SM22]*mDet2_34_01;
0300   const Scalar mDet3_234_013 = rhs[SM20]*mDet2_34_13 - rhs[SM21]*mDet2_34_03 + rhs[SM23]*mDet2_34_01;
0301   const Scalar mDet3_234_014 = rhs[SM20]*mDet2_34_14 - rhs[SM21]*mDet2_34_04 + rhs[SM24]*mDet2_34_01;
0302   const Scalar mDet3_234_023 = rhs[SM20]*mDet2_34_23 - rhs[SM22]*mDet2_34_03 + rhs[SM23]*mDet2_34_02;
0303   const Scalar mDet3_234_024 = rhs[SM20]*mDet2_34_24 - rhs[SM22]*mDet2_34_04 + rhs[SM24]*mDet2_34_02;
0304   const Scalar mDet3_234_034 = rhs[SM20]*mDet2_34_34 - rhs[SM23]*mDet2_34_04 + rhs[SM24]*mDet2_34_03;
0305   const Scalar mDet3_234_123 = rhs[SM21]*mDet2_34_23 - rhs[SM22]*mDet2_34_13 + rhs[SM23]*mDet2_34_12;
0306   const Scalar mDet3_234_124 = rhs[SM21]*mDet2_34_24 - rhs[SM22]*mDet2_34_14 + rhs[SM24]*mDet2_34_12;
0307   const Scalar mDet3_234_134 = rhs[SM21]*mDet2_34_34 - rhs[SM23]*mDet2_34_14 + rhs[SM24]*mDet2_34_13;
0308   const Scalar mDet3_234_234 = rhs[SM22]*mDet2_34_34 - rhs[SM23]*mDet2_34_24 + rhs[SM24]*mDet2_34_23;
0309
0310   // Find all NECESSARY 4x4 dets:   (15 of them)
0311
0312   const Scalar mDet4_0123_0123 = rhs[SM00]*mDet3_123_123 - rhs[SM01]*mDet3_123_023
0313                                  + rhs[SM02]*mDet3_123_013 - rhs[SM03]*mDet3_123_012;
0314   const Scalar mDet4_0124_0123 = rhs[SM00]*mDet3_124_123 - rhs[SM01]*mDet3_124_023
0315                                  + rhs[SM02]*mDet3_124_013 - rhs[SM03]*mDet3_124_012;
0316   const Scalar mDet4_0124_0124 = rhs[SM00]*mDet3_124_124 - rhs[SM01]*mDet3_124_024
0317                                  + rhs[SM02]*mDet3_124_014 - rhs[SM04]*mDet3_124_012;
0318   const Scalar mDet4_0134_0123 = rhs[SM00]*mDet3_134_123 - rhs[SM01]*mDet3_134_023
0319                                  + rhs[SM02]*mDet3_134_013 - rhs[SM03]*mDet3_134_012;
0320   const Scalar mDet4_0134_0124 = rhs[SM00]*mDet3_134_124 - rhs[SM01]*mDet3_134_024
0321                                  + rhs[SM02]*mDet3_134_014 - rhs[SM04]*mDet3_134_012;
0322   const Scalar mDet4_0134_0134 = rhs[SM00]*mDet3_134_134 - rhs[SM01]*mDet3_134_034
0323                                  + rhs[SM03]*mDet3_134_014 - rhs[SM04]*mDet3_134_013;
0324   const Scalar mDet4_0234_0123 = rhs[SM00]*mDet3_234_123 - rhs[SM01]*mDet3_234_023
0325                                  + rhs[SM02]*mDet3_234_013 - rhs[SM03]*mDet3_234_012;
0326   const Scalar mDet4_0234_0124 = rhs[SM00]*mDet3_234_124 - rhs[SM01]*mDet3_234_024
0327                                  + rhs[SM02]*mDet3_234_014 - rhs[SM04]*mDet3_234_012;
0328   const Scalar mDet4_0234_0134 = rhs[SM00]*mDet3_234_134 - rhs[SM01]*mDet3_234_034
0329                                  + rhs[SM03]*mDet3_234_014 - rhs[SM04]*mDet3_234_013;
0330   const Scalar mDet4_0234_0234 = rhs[SM00]*mDet3_234_234 - rhs[SM02]*mDet3_234_034
0331                                  + rhs[SM03]*mDet3_234_024 - rhs[SM04]*mDet3_234_023;
0332   const Scalar mDet4_1234_0123 = rhs[SM10]*mDet3_234_123 - rhs[SM11]*mDet3_234_023
0333                                  + rhs[SM12]*mDet3_234_013 - rhs[SM13]*mDet3_234_012;
0334   const Scalar mDet4_1234_0124 = rhs[SM10]*mDet3_234_124 - rhs[SM11]*mDet3_234_024
0335                                  + rhs[SM12]*mDet3_234_014 - rhs[SM14]*mDet3_234_012;
0336   const Scalar mDet4_1234_0134 = rhs[SM10]*mDet3_234_134 - rhs[SM11]*mDet3_234_034
0337                                  + rhs[SM13]*mDet3_234_014 - rhs[SM14]*mDet3_234_013;
0338   const Scalar mDet4_1234_0234 = rhs[SM10]*mDet3_234_234 - rhs[SM12]*mDet3_234_034
0339                                  + rhs[SM13]*mDet3_234_024 - rhs[SM14]*mDet3_234_023;
0340   const Scalar mDet4_1234_1234 = rhs[SM11]*mDet3_234_234 - rhs[SM12]*mDet3_234_134
0341                                  + rhs[SM13]*mDet3_234_124 - rhs[SM14]*mDet3_234_123;
0342
0343   // Find the 5x5 det:
0344
0345   const Scalar det = rhs[SM00]*mDet4_1234_1234 - rhs[SM01]*mDet4_1234_0234 + rhs[SM02]*mDet4_1234_0134
0346                        - rhs[SM03]*mDet4_1234_0124 + rhs[SM04]*mDet4_1234_0123;
0347 //   if (determ)
0348 //     *determ = det;
0349
0350   if ( det == 0 )
0351     return false;
0352
0353   const Scalar oneOverDet = 1.0f / det;
0354   const Scalar mn1OverDet = - oneOverDet;
0355
0356   rhs[SM00] = mDet4_1234_1234 * oneOverDet;
0357   rhs[SM01] = mDet4_1234_0234 * mn1OverDet;
0358   rhs[SM02] = mDet4_1234_0134 * oneOverDet;
0359   rhs[SM03] = mDet4_1234_0124 * mn1OverDet;
0360   rhs[SM04] = mDet4_1234_0123 * oneOverDet;
0361
0362   rhs[SM11] = mDet4_0234_0234 * oneOverDet;
0363   rhs[SM12] = mDet4_0234_0134 * mn1OverDet;
0364   rhs[SM13] = mDet4_0234_0124 * oneOverDet;
0365   rhs[SM14] = mDet4_0234_0123 * mn1OverDet;
0366
0367   rhs[SM22] = mDet4_0134_0134 * oneOverDet;
0368   rhs[SM23] = mDet4_0134_0124 * mn1OverDet;
0369   rhs[SM24] = mDet4_0134_0123 * oneOverDet;
0370
0371   rhs[SM33] = mDet4_0124_0124 * oneOverDet;
0372   rhs[SM34] = mDet4_0124_0123 * mn1OverDet;
0373
0374   rhs[SM44] = mDet4_0123_0123 * oneOverDet;
0375
0376
0377   return true;
0378 }
0379
0380
0381
0382   }  // namespace Math
0383
0384 }  // namespace ROOT
0385
0386
0387 // undef Smacros to avoid conflicts
0388
0389 // undef SF 4x4
0390 #undef SF00
0391 #undef SF01
0392 #undef SF02
0393 #undef SF03
0394
0395 #undef SF10
0396 #undef SF11
0397 #undef SF12
0398 #undef SF13
0399
0400 #undef SF20
0401 #undef SF21
0402 #undef SF22
0403 #undef SF23
0404
0405 #undef SF30
0406 #undef SF31
0407 #undef SF32
0408 #undef SF33
0409
0410 // undef SM 5x5
0411 #undef SM00
0412 #undef SM01
0413 #undef SM02
0414 #undef SM03
0415 #undef SM04
0416
0417 #undef SM10
0418 #undef SM11
0419 #undef SM12
0420 #undef SM13
0421 #undef SM14
0422
0423 #undef SM20
0424 #undef SM21
0425 #undef SM22
0426 #undef SM23
0427 #undef SM24
0428
0429 #undef SM30
0430 #undef SM31
0431 #undef SM32
0432 #undef SM33
0433 #undef SM34
0434
0435 #undef SM40
0436 #undef SM41
0437 #undef SM42
0438 #undef SM43
0439 #undef SM44
0440
0441
0442
0443
0444 #endif