Why?

ROOT
1

why do I get 2 different results for calling a constructor looking the same way to me. I can steer it with placing brackets, but I really would like to understand what the real difference between

Transform_global_to_lmd_local_ref and Transform_global_to_lmd_local_val is.

Could you please check the simple example I have attached?

Since it is so simple I place it also inline for others readers:

// to compile 
// g++ $(bash root-config --cflags --glibs) -o test_prog test_prog.cxx
// then run 
// ./test_prog
// why do I get 2 different results?

#include "TMatrixD.h"

TMatrixD Get_rotation_matrix(){
	TMatrixD result(3,3);
	result[0][0] = 2.;
	result[1][0] = -1.5;
	result[2][0] = 3.;
	result[0][1] = 6.;
	result[1][1] = -2.5;
	result[2][1] = 2.5;
	result[0][2] = 10.;
	result[1][2] = -0.5;
	result[2][2] = 1.;
	return result;
}

TMatrixD Get_matrix_for_rotation(){
	TMatrixD result(3,3);
	result[0][0] = 1.;
	result[1][0] = 2.5;
	result[2][0] = -2.;
	result[0][1] = -1.;
	result[1][1] = -3.5;
	result[2][1] = 0.5;
	result[0][2] = 4.;
	result[1][2] = 0.5;
	result[2][2] = -4.;
	return result;
}

TMatrixD& Transform_global_to_lmd_local_ref(const TMatrixD& matrix){
	TMatrixD rotmatrix(Get_rotation_matrix());
	return *(new TMatrixD(rotmatrix*matrix*rotmatrix.T()));
}

TMatrixD Transform_global_to_lmd_local_val(const TMatrixD& matrix){
	TMatrixD rotmatrix(Get_rotation_matrix());
	return TMatrixD(rotmatrix*matrix*rotmatrix.T());
}

TMatrixD Transform_global_to_lmd_local_ana(const TMatrixD& matrix){
	TMatrixD rotmatrix(Get_rotation_matrix());
	return TMatrixD(rotmatrix*TMatrixD(matrix,TMatrixD::kMultTranspose,rotmatrix));
}

int main(void)
{
	TMatrixD matrix(Get_matrix_for_rotation());
	TMatrixD& rotated_matrix1 = Transform_global_to_lmd_local_ref(matrix);
	const TMatrixD& rotated_matrix2 = Transform_global_to_lmd_local_val(matrix);
	const TMatrixD& rotated_matrix3 = Transform_global_to_lmd_local_ana(matrix);

	rotated_matrix1.Print();
	rotated_matrix2.Print();
	rotated_matrix3.Print();

	return 0;
}

In my case I get

[jasinski@PALUMA-121121 temp]$ ./test_prog 

3x3 matrix is as follows

     |      0    |      1    |      2    |
--------------------------------------------
   0 |       -404          64         -83 
   1 |         -6      -12.25        -0.5 
   2 |       39.5       12.62       2.875 


3x3 matrix is as follows

     |      0    |      1    |      2    |
--------------------------------------------
   0 |      -34.5       2.125      -16.12 
   1 |        141        -8.5           7 
   2 |      324.5      -8.625       36.62 


3x3 matrix is as follows

     |      0    |      1    |      2    |
--------------------------------------------
   0 |       -404          64         -83 
   1 |         -6      -12.25        -0.5 
   2 |       39.5       12.62       2.875

as a result and I really don’t understand why the second one differs from the first one.
test_prog.cxx (1.58 KB)

2

I just checked, there is no way to steer Transform_global_to_lmd_local_val with brackets, so it give the same result as Transform_global_to_lmd_local_ref.

return TMatrixD(rotmatrix*(matrix*rotmatrix.T()));

and

return TMatrixD((rotmatrix*matrix)*rotmatrix.T());

give the same result.

for further information:
root is 5.34/01 and gcc is 4.7.2

Sorry for being so stupid that I do not understand this

3

Ok, still we don’t understand why the order of executions is different but at least we realized that

TMatrixD::T()

actually changes the entries of the matrix. Since it is a function returning the matrix I thought it would be a copy of the matrix. Is there a way to rotate a matrix in an performant and elegant looking way?

4

The matrix package is supposed to not allocate more memory without
explicit user action, so the matrix is transposed in place.

Obviously the order of execution is different in both cases and will
have different final results ( because of the in place transpose).

In stead of the allocating a new matrix and transposing it,
you will have to do it:

(TMatrixD(rotmatrix)).T()

TMatrixD(rotmatrix,TMatrixD::kTransposed)

btw. Why are you wrapping the return objects in a TMatrixD() statement ?