removed c0-1 multiplication.

[phys/u1casc.git] / u1casc-ordinary / obs_diagcorr.hpp

diff --git a/u1casc-ordinary/obs_diagcorr.hpp b/u1casc-ordinary/obs_diagcorr.hpp
index f51df7cd31fb23b7d9af3ff826911f66c6368ebb..6666b6c5938c9f73a18a0a294db19083a66a3f13 100644 (file)
--- a/u1casc-ordinary/obs_diagcorr.hpp
+++ b/u1casc-ordinary/obs_diagcorr.hpp
@@ -47,8 +47,6 @@ private:
   gsl_vector_complex **measurements_disco;
   
   static void cdiag (gsl_vector *v, gsl_matrix_complex *m);
-
-  static void invertcm (gsl_matrix_complex *invm, const gsl_matrix_complex *m);
 };
 
 obs_diagcorr::obs_diagcorr(o815 *_O815) : o815::obs("diagcorr", 
@@ -97,20 +95,6 @@ void obs_diagcorr::_meas(bool loadedobs, const int& nthmeas) {
   }
 };
 
-void obs_diagcorr::invertcm (gsl_matrix_complex *invm, const gsl_matrix_complex *m)
-{
-  int s;
-  gsl_permutation * perm = gsl_permutation_alloc (m->size1);
-  
-  gsl_matrix_complex *mtmp = gsl_matrix_complex_alloc(m->size1,m->size2);
-  gsl_matrix_complex_memcpy(mtmp, m);
-
-  gsl_linalg_complex_LU_decomp (mtmp, perm, &s);
-  gsl_linalg_complex_LU_invert (mtmp, perm, invm);
-
-  gsl_matrix_complex_free(mtmp);
-}
-
 void obs_diagcorr::cdiag (gsl_vector *v, gsl_matrix_complex *m)
 {
   gsl_matrix_complex *evec = gsl_matrix_complex_alloc(4,4);  
@@ -126,14 +110,12 @@ void obs_diagcorr::cdiag (gsl_vector *v, gsl_matrix_complex *m)
 void obs_diagcorr::_finish() {
   gsl_matrix_complex *totalval = gsl_matrix_complex_alloc(4,4);
   gsl_matrix_complex *tmpmatrix = gsl_matrix_complex_alloc(4,4);
-  gsl_matrix_complex *tmpmatrix2 = gsl_matrix_complex_alloc(4,4);
   gsl_vector *tmpvec = gsl_vector_alloc(4);
   gsl_vector_complex *tmpvecc = gsl_vector_complex_alloc(4);
   gsl_vector *tmpvec2 = gsl_vector_alloc(4);
   gsl_vector *jackres = gsl_vector_alloc(4);
   gsl_vector *jackerrornorm = gsl_vector_alloc(4);
   gsl_vector_complex *totaldisco = gsl_vector_complex_alloc(4);
-  gsl_matrix_complex *invC0 = NULL;
   
   gsl_vector_complex_set_zero(totaldisco);
   for (int imeas=0; imeas<O815->comargs.nmeas; imeas++)
@@ -160,15 +142,7 @@ void obs_diagcorr::_finish() {
                                                 )
                                );
       }
-
-    if (invC0 == NULL) {
-      invC0 = gsl_matrix_complex_alloc(4,4);
-      invertcm(invC0, tmpmatrix); 
-    }
-
-    gsl_blas_zgemm (CblasNoTrans, CblasNoTrans, gsl_complex_rect(1.0,0.0),
-                   invC0, tmpmatrix, gsl_complex_rect(0.0,0.0), tmpmatrix2);
-    cdiag(jackres, tmpmatrix2);
+    cdiag(jackres, tmpmatrix);
 
     for (int imeas=0; imeas<O815->comargs.nmeas; imeas++) {
       gsl_matrix_complex_memcpy (tmpmatrix, totalval);
@@ -187,12 +161,8 @@ void obs_diagcorr::_finish() {
                                                   discopart
                                                   )
                                  );
-       }
-
-      gsl_blas_zgemm (CblasNoTrans, CblasNoTrans, gsl_complex_rect(1.0,0.0),
-                     invC0, tmpmatrix, gsl_complex_rect(0.0,0.0), tmpmatrix2);
-      
-      cdiag(tmpvec, tmpmatrix2);
+       }      
+      cdiag(tmpvec, tmpmatrix);
       gsl_vector_sub(tmpvec, jackres);
       gsl_vector_memcpy(tmpvec2, tmpvec);
       gsl_vector_mul(tmpvec, tmpvec2);