X-Git-Url: https://pintos-os.org/cgi-bin/gitweb.cgi?a=blobdiff_plain;f=src%2Fmath%2Fcovariance.c;h=79fce25cef28362dbe7b404377d436e92727111c;hb=bf868380ccff985a9eeb4e40f307118548411e7c;hp=a8c71dc0b357289a65751814f0e5c9499ac73478;hpb=6e097c89af440da90b43ce90864394c4d0c843d5;p=pspp

diff --git a/src/math/covariance.c b/src/math/covariance.c
index a8c71dc0b3..79fce25cef 100644
--- a/src/math/covariance.c
+++ b/src/math/covariance.c
@@ -70,6 +70,9 @@ resize_matrix (gsl_matrix *in, size_t new_size)
 
 struct covariance
 {
+  /* True if the covariances are centerered. (ie Real covariances) */
+  bool centered;
+
   /* The variables for which the covariance matrix is to be calculated. */
   size_t n_vars;
   const struct variable *const *vars;
@@ -126,7 +129,7 @@ struct covariance
    be identical.  If missing values are involved, then element (i,j)
    is the moment of the i th variable, when paired with the j th variable.
  */
-const gsl_matrix *
+gsl_matrix *
 covariance_moments (const struct covariance *cov, int m)
 {
   return cov->moments[m];
@@ -138,11 +141,13 @@ covariance_moments (const struct covariance *cov, int m)
  */
 struct covariance *
 covariance_1pass_create (size_t n_vars, const struct variable *const *vars,
-			 const struct variable *weight, enum mv_class exclude)
+			 const struct variable *weight, enum mv_class exclude,
+			 bool centered)
 {
   size_t i;
-  struct covariance *cov = xzalloc (sizeof *cov);
+  struct covariance *cov = XZALLOC (struct covariance);
 
+  cov->centered = centered;
   cov->passes = 1;
   cov->state = 0;
   cov->pass_one_first_case_seen = cov->pass_two_first_case_seen = false;
@@ -160,7 +165,7 @@ covariance_1pass_create (size_t n_vars, const struct variable *const *vars,
 
   cov->exclude = exclude;
 
-  cov->n_cm = (n_vars * (n_vars - 1)  ) / 2;
+  cov->n_cm = (n_vars * (n_vars - 1)) / 2;
 
 
   cov->cm = xcalloc (cov->n_cm, sizeof *cov->cm);
@@ -178,11 +183,13 @@ covariance_1pass_create (size_t n_vars, const struct variable *const *vars,
 struct covariance *
 covariance_2pass_create (size_t n_vars, const struct variable *const *vars,
 			 struct categoricals *cats,
-			 const struct variable *wv, enum mv_class exclude)
+			 const struct variable *wv, enum mv_class exclude,
+			 bool centered)
 {
   size_t i;
   struct covariance *cov = xmalloc (sizeof *cov);
 
+  cov->centered = centered;
   cov->passes = 2;
   cov->state = 0;
   cov->pass_one_first_case_seen = cov->pass_two_first_case_seen = false;
@@ -225,13 +232,13 @@ cm_idx (const struct covariance *cov, int i, int j)
   assert (i >= 0);
   assert (j < cov->dim);
 
-  if ( i == 0)
+  if (i == 0)
     return -1;
 
   if (j >= cov->dim - 1)
     return -1;
 
-  if ( i <= j)
+  if (i <= j)
     return -1 ;
 
   as = nj * (nj + 1) ;
@@ -255,14 +262,14 @@ is_missing (const struct covariance *cov, int i, const struct ccase *c)
 
   const union value *val = case_data (c, var);
 
-  return var_is_value_missing (var, val, cov->exclude);
+  return (var_is_value_missing (var, val) & cov->exclude) != 0;
 }
 
 
 static double
 get_val (const struct covariance *cov, int i, const struct ccase *c)
 {
-  if ( i < cov->n_vars)
+  if (i < cov->n_vars)
     {
       const struct variable *var = cov->vars[i];
 
@@ -294,7 +301,7 @@ void
 covariance_accumulate_pass1 (struct covariance *cov, const struct ccase *c)
 {
   size_t i, j, m;
-  const double weight = cov->wv ? case_data (c, cov->wv)->f : 1.0;
+  const double weight = cov->wv ? case_num (c, cov->wv) : 1.0;
 
   assert (cov->passes == 2);
   if (!cov->pass_one_first_case_seen)
@@ -310,14 +317,14 @@ covariance_accumulate_pass1 (struct covariance *cov, const struct ccase *c)
     {
       double v1 = get_val (cov, i, c);
 
-      if ( is_missing (cov, i, c))
+      if (is_missing (cov, i, c))
 	continue;
 
       for (j = 0 ; j < cov->dim; ++j)
 	{
 	  double pwr = 1.0;
 
-	  if ( is_missing (cov, j, c))
+	  if (is_missing (cov, j, c))
 	    continue;
 
 	  for (m = 0 ; m <= MOMENT_MEAN; ++m)
@@ -339,7 +346,7 @@ void
 covariance_accumulate_pass2 (struct covariance *cov, const struct ccase *c)
 {
   size_t i, j;
-  const double weight = cov->wv ? case_data (c, cov->wv)->f : 1.0;
+  const double weight = cov->wv ? case_num (c, cov->wv) : 1.0;
 
   assert (cov->passes == 2);
   assert (cov->state >= 1);
@@ -358,7 +365,7 @@ covariance_accumulate_pass2 (struct covariance *cov, const struct ccase *c)
       if (cov->categoricals)
 	cov->dim += categoricals_df_total (cov->categoricals);
 
-      cov->n_cm = (cov->dim * (cov->dim - 1)  ) / 2;
+      cov->n_cm = (cov->dim * (cov->dim - 1)) / 2;
       cov->cm = xcalloc (cov->n_cm, sizeof *cov->cm);
 
       /* Grow the moment matrices so that they're large enough to accommodate the
@@ -412,7 +419,7 @@ covariance_accumulate_pass2 (struct covariance *cov, const struct ccase *c)
     {
       double v1 = get_val (cov, i, c);
 
-      if ( is_missing (cov, i, c))
+      if (is_missing (cov, i, c))
 	continue;
 
       for (j = 0 ; j < cov->dim; ++j)
@@ -423,7 +430,7 @@ covariance_accumulate_pass2 (struct covariance *cov, const struct ccase *c)
 
 	  const double s = pow2 (v1 - gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)) * weight;
 
-	  if ( is_missing (cov, j, c))
+	  if (is_missing (cov, j, c))
 	    continue;
 
 	  {
@@ -457,13 +464,13 @@ void
 covariance_accumulate (struct covariance *cov, const struct ccase *c)
 {
   size_t i, j, m;
-  const double weight = cov->wv ? case_data (c, cov->wv)->f : 1.0;
+  const double weight = cov->wv ? case_num (c, cov->wv) : 1.0;
 
   assert (cov->passes == 1);
 
-  if ( !cov->pass_one_first_case_seen)
+  if (!cov->pass_one_first_case_seen)
     {
-      assert ( cov->state == 0);
+      assert (cov->state == 0);
       cov->state = 1;
     }
 
@@ -471,7 +478,7 @@ covariance_accumulate (struct covariance *cov, const struct ccase *c)
     {
       const union value *val1 = case_data (c, cov->vars[i]);
 
-      if ( is_missing (cov, i, c))
+      if (is_missing (cov, i, c))
 	continue;
 
       for (j = 0 ; j < cov->dim; ++j)
@@ -480,7 +487,7 @@ covariance_accumulate (struct covariance *cov, const struct ccase *c)
 	  int idx;
 	  const union value *val2 = case_data (c, cov->vars[j]);
 
-	  if ( is_missing (cov, j, c))
+	  if (is_missing (cov, j, c))
 	    continue;
 
 	  idx = cm_idx (cov, i, j);
@@ -514,7 +521,7 @@ cm_to_gsl (struct covariance *cov)
   gsl_matrix *m = gsl_matrix_calloc (cov->dim, cov->dim);
 
   /* Copy the non-diagonal elements from cov->cm */
-  for ( j = 0 ; j < cov->dim - 1; ++j)
+  for (j = 0 ; j < cov->dim - 1; ++j)
     {
       for (i = j+1 ; i < cov->dim; ++i)
 	{
@@ -548,7 +555,7 @@ covariance_calculate_double_pass (struct covariance *cov)
 	  *x /= gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
 
 	  idx = cm_idx (cov, i, j);
-	  if ( idx >= 0)
+	  if (idx >= 0)
 	    {
 	      x = &cov->cm [idx];
 	      *x /= gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
@@ -568,7 +575,7 @@ covariance_calculate_single_pass (struct covariance *cov)
   for (m = 0; m < n_MOMENTS; ++m)
     {
       /* Divide the moments by the number of samples */
-      if ( m > 0)
+      if (m > 0)
 	{
 	  for (i = 0 ; i < cov->dim; ++i)
 	    {
@@ -577,26 +584,29 @@ covariance_calculate_single_pass (struct covariance *cov)
 		  double *x = gsl_matrix_ptr (cov->moments[m], i, j);
 		  *x /= gsl_matrix_get (cov->moments[0], i, j);
 
-		  if ( m == MOMENT_VARIANCE)
+		  if (m == MOMENT_VARIANCE)
 		    *x -= pow2 (gsl_matrix_get (cov->moments[1], i, j));
 		}
 	    }
 	}
     }
 
-  /* Centre the moments */
-  for ( j = 0 ; j < cov->dim - 1; ++j)
+  if (cov->centered)
     {
-      for (i = j + 1 ; i < cov->dim; ++i)
+      /* Centre the moments */
+      for (j = 0 ; j < cov->dim - 1; ++j)
 	{
-	  double *x = &cov->cm [cm_idx (cov, i, j)];
+	  for (i = j + 1 ; i < cov->dim; ++i)
+	    {
+	      double *x = &cov->cm [cm_idx (cov, i, j)];
 
-	  *x /= gsl_matrix_get (cov->moments[0], i, j);
+	      *x /= gsl_matrix_get (cov->moments[0], i, j);
 
-	  *x -=
-	    gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
-	    *
-	    gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i);
+	      *x -=
+		gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
+		*
+		gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i);
+	    }
 	}
     }
 
@@ -612,7 +622,7 @@ covariance_calculate_single_pass (struct covariance *cov)
 gsl_matrix *
 covariance_calculate (struct covariance *cov)
 {
-  if ( cov->state <= 0 )
+  if (cov->state <= 0)
     return NULL;
 
   switch (cov->passes)
@@ -642,26 +652,30 @@ covariance_calculate_single_pass_unnormalized (struct covariance *cov)
 {
   size_t i, j;
 
-  for (i = 0 ; i < cov->dim; ++i)
+  if (cov->centered)
     {
-      for (j = 0 ; j < cov->dim; ++j)
+      for (i = 0 ; i < cov->dim; ++i)
 	{
-	  double *x = gsl_matrix_ptr (cov->moments[MOMENT_VARIANCE], i, j);
-	  *x -= pow2 (gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j))
-	    / gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+	  for (j = 0 ; j < cov->dim; ++j)
+	    {
+	      double *x = gsl_matrix_ptr (cov->moments[MOMENT_VARIANCE], i, j);
+	      *x -= pow2 (gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j))
+		/ gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+	    }
 	}
-    }
-  for ( j = 0 ; j < cov->dim - 1; ++j)
-    {
-      for (i = j + 1 ; i < cov->dim; ++i)
+
+      for (j = 0 ; j < cov->dim - 1; ++j)
 	{
-	  double *x = &cov->cm [cm_idx (cov, i, j)];
+	  for (i = j + 1 ; i < cov->dim; ++i)
+	    {
+	      double *x = &cov->cm [cm_idx (cov, i, j)];
 
-	  *x -=
-	    gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
-	    *
-	    gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i)
-	  / gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+	      *x -=
+		gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
+		*
+		gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i)
+		/ gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+	    }
 	}
     }
 
@@ -676,7 +690,7 @@ covariance_calculate_single_pass_unnormalized (struct covariance *cov)
 const gsl_matrix *
 covariance_calculate_unnormalized (struct covariance *cov)
 {
-  if ( cov->state <= 0 )
+  if (cov->state <= 0)
     return NULL;
 
   if (cov->unnormalised != NULL)
@@ -740,62 +754,46 @@ covariance_dim (const struct covariance * cov)
 */
 
 #include "libpspp/str.h"
-#include "output/tab.h"
+#include "output/pivot-table.h"
 #include "data/format.h"
 
 
 /* Create a table which can be populated with the encodings for
    the covariance matrix COV */
-struct tab_table *
-covariance_dump_enc_header (const struct covariance *cov, int length)
+struct pivot_table *
+covariance_dump_enc_header (const struct covariance *cov)
 {
-  struct tab_table *t = tab_create (cov->dim,  length);
-  int n;
-  int i;
-
-  tab_title (t, "Covariance Encoding");
-
-  tab_box (t,
-	   TAL_2, TAL_2, 0, 0,
-	   0, 0,   tab_nc (t) - 1,   tab_nr (t) - 1);
-
-  tab_hline (t, TAL_2, 0, tab_nc (t) - 1, 1);
-
-
-  for (i = 0 ; i < cov->n_vars; ++i)
+  struct pivot_table *table = pivot_table_create ("Covariance Encoding");
+
+  struct pivot_dimension *factors = pivot_dimension_create (
+    table, PIVOT_AXIS_COLUMN, "Factor");
+  for (size_t i = 0 ; i < cov->n_vars; ++i)
+    pivot_category_create_leaf (factors->root,
+                                pivot_value_new_variable (cov->vars[i]));
+  for (size_t i = 0, n = 0; i < cov->dim - cov->n_vars; n++)
     {
-      tab_text (t, i, 0, TAT_TITLE, var_get_name (cov->vars[i]));
-      tab_vline (t, TAL_1, i + 1, 0, tab_nr (t) - 1);
-    }
-
-  n = 0;
-  while (i < cov->dim)
-    {
-      struct string str;
-      int idx = i - cov->n_vars;
       const struct interaction *iact =
-	categoricals_get_interaction_by_subscript (cov->categoricals, idx);
-      int df;
+	categoricals_get_interaction_by_subscript (cov->categoricals, i);
 
-      ds_init_empty (&str);
+      struct string str = DS_EMPTY_INITIALIZER;
       interaction_to_string (iact, &str);
+      struct pivot_category *group = pivot_category_create_group__ (
+        factors->root,
+        pivot_value_new_user_text_nocopy (ds_steal_cstr (&str)));
 
-      df = categoricals_df (cov->categoricals, n);
-
-      tab_joint_text (t,
-		      i, 0,
-		      i + df - 1, 0,
-		      TAT_TITLE, ds_cstr (&str));
-
-      if (i + df < tab_nr (t) - 1)
-	tab_vline (t, TAL_1, i + df, 0, tab_nr (t) - 1);
+      int df = categoricals_df (cov->categoricals, n);
+      for (int j = 0; j < df; j++)
+        pivot_category_create_leaf_rc (group, pivot_value_new_integer (j),
+                                       PIVOT_RC_INTEGER);
 
       i += df;
-      n++;
-      ds_destroy (&str);
     }
 
-  return t;
+  struct pivot_dimension *matrix = pivot_dimension_create (
+    table, PIVOT_AXIS_ROW, "Matrix", "Matrix");
+  matrix->hide_all_labels = true;
+
+  return table;
 }
 
 
@@ -805,14 +803,13 @@ covariance_dump_enc_header (const struct covariance *cov, int length)
  */
 void
 covariance_dump_enc (const struct covariance *cov, const struct ccase *c,
-		     struct tab_table *t)
+		     struct pivot_table *table)
 {
-  static int row = 0;
-  int i;
-  ++row;
-  for (i = 0 ; i < cov->dim; ++i)
-    {
-      double v = get_val (cov, i, c);
-      tab_double (t, i, row, 0, v, i < cov->n_vars ? NULL : &F_8_0, RC_OTHER);
-    }
+  int row = pivot_category_create_leaf (
+    table->dimensions[1]->root,
+    pivot_value_new_integer (table->dimensions[1]->n_leaves));
+
+  for (int i = 0 ; i < cov->dim; ++i)
+    pivot_table_put2 (
+      table, i, row, pivot_value_new_number (get_val (cov, i, c)));
 }