wrote old auton-style multi-tree n-point, did very basic tests, needs more testing.

2011-02-09 21:45:56 +00:00
parent 915cd9e4fc
commit 5adef387dc
5 changed files with 306 additions and 4 deletions
@@ -9,6 +9,8 @@ set(SOURCES
 	permutations.cc
 	single_matcher.h
 	single_matcher.cc
+	single_bandwidth_alg.h
+	single_bandwidth_alg.cc
 )

 # add directory name to sources
@@ -7,8 +7,7 @@
 *
 */

-// take this out later, just debugging with it
-#include "single_matcher.h"
+#include "single_bandwidth_alg.h"

 using namespace npt;

@@ -37,12 +36,27 @@ int main(int argc, char* argv[]) {
  
  upper_bds.load("test_upper_bds.csv");
  
+  arma::mat data;
+  data.load("test_data.csv");
+  
+  arma::colvec weights(data.n_cols);
+  weights.fill(1.0);
+  
+  index_t leaf_size = 1;
+  
+  SingleBandwidthAlg alg(data, weights, data.n_cols, leaf_size, lower_bds, 
+                         upper_bds);
+  
+  std::cout << "Num tuples: " << alg.num_tuples() << "\n";
+  
+  /*
  SingleMatcher matcher(3, lower_bds, upper_bds);
  
  std::vector<bool> perm_ok1(6, true);
  std::cout << "Testing points: " << matcher.TestPointPair(0.5, 0, 1, perm_ok1)
  << "\n";
-
+  */
+  
  //std::vector<bool> perm_ok2(6, true);
  //std::cout << "Testing boxes: " << TestHrectPair(5.0, );
  
@@ -9,3 +9,216 @@

 #include "single_bandwidth_alg.h"

+bool npt::SingleBandwidthAlg::CheckNodeList_(std::vector<SingleNode*>& nodes) {
+  
+  bool can_prune = false;
+  
+  // note that this has to enforce symmetry
+  std::vector<bool> permutation_ok(matcher_.num_permutations(), true);
+  
+  // iterate over all nodes
+  // IMPORTANT: right now, I'm exiting when I can prune
+  // I need to double check that this works
+  for (index_t i = 0; !can_prune && i < tuple_size_; i++) {
+    
+    SingleNode* node_i = nodes[i];
+    
+    // iterate over all nodes > i
+    for (index_t j = i+1; !can_prune && j < tuple_size_; j++) {
+      
+      SingleNode* node_j = nodes[j];
+      
+      // check for symmetry
+      if (node_j->end() <= node_i->begin()) {
+        return false;
+      } // symmetry check
+      
+      can_prune = matcher_.TestHrectPair(node_i->bound(), node_j->bound(),
+                                         i, j, permutation_ok);
+      
+    } // for j
+    
+  } // for i
+  
+  return can_prune;
+  
+} // CheckNodeList
+
+
+
+void npt::SingleBandwidthAlg::BaseCaseHelper_(std::vector<std::vector<index_t> >& point_sets,
+                                              std::vector<bool>& permutation_ok,
+                                              std::vector<index_t>& points_in_tuple,
+                                              int k) {
+  
+
+  std::vector<bool> permutation_ok_copy(permutation_ok);
+
+  bool bad_symmetry = false;
+  
+  std::vector<index_t>& k_rows = point_sets[k];
+  
+  // iterate over possible kth members of the tuple
+  for (index_t i = 0; i < k_rows.size(); i++) {
+    
+    index_t point_i_index = k_rows[i];
+    bool this_point_works = true;
+    
+    arma::colvec vec_i = data_points_.col(point_i_index);
+    
+    // TODO: Does this leak memory?
+    permutation_ok_copy.assign(permutation_ok.begin(), permutation_ok.end());
+    
+    // loop over points already in the tuple and check against them
+    for (index_t j = 0; this_point_works && j < k; j++) {
+      
+      index_t point_j_index = points_in_tuple[j];
+
+      // j comes before i in the tuple, so it should have a lower index
+      bad_symmetry = (point_i_index <= point_j_index);
+      
+      if (!bad_symmetry) {
+        
+        arma::colvec vec_j = data_points_.col(j);
+        
+        double point_dist_sq = la::DistanceSqEuclidean(vec_i, vec_j);
+        
+        this_point_works = matcher_.TestPointPair(point_dist_sq, j, k, 
+                                                  permutation_ok_copy);
+        
+      } // check the distances across permutations
+      
+    } // for j
+    
+    // point i fits in the tuple
+    if (this_point_works && !bad_symmetry) {
+      
+      points_in_tuple[k] = point_i_index;
+      
+      // are we finished?
+      if (k == tuple_size_ - 1) {
+        
+        num_tuples_++;
+        
+        double this_weight = 1.0;
+        
+        for (index_t tuple_ind = 0; tuple_ind < tuple_size_; tuple_ind++) {
+          
+          this_weight *= data_weights_(points_in_tuple[tuple_ind]);
+          
+        } // iterate over the tuple
+        
+        weighted_num_tuples_ += this_weight;
+        
+      } 
+      else {
+        
+        BaseCaseHelper_(point_sets, permutation_ok_copy, points_in_tuple, k+1);
+        
+      } // need to add more points to finish the tuple
+      
+    } // point i fits
+    
+  } // for i
+  
+} // BaseCaseHelper_
+
+
+
+void npt::SingleBandwidthAlg::BaseCase_(std::vector<SingleNode*>& nodes) {
+  
+  std::vector<std::vector<index_t> > point_sets;
+
+  // TODO: can this be done more efficiently?
+  
+  // Make a 2D array of the points in the nodes 
+  // iterate over nodes
+  for (index_t node_ind = 0; node_ind < tuple_size_; node_ind++) {
+    
+    point_sets[node_ind].resize(nodes[node_ind]->count());
+    
+    // fill in points in the node
+    for (index_t point_ind = nodes[node_ind]->begin(); 
+         point_ind < nodes[node_ind]->end(); point_ind++) {
+      
+      point_sets[node_ind][point_ind] = point_ind;
+      
+    } // points
+    
+  } // nodes
+  
+  std::vector<bool> permutation_ok(matcher_.num_permutations(), true);
+  
+  std::vector<index_t> points_in_tuple(tuple_size_, -1);
+  
+  BaseCaseHelper_(point_sets, permutation_ok, points_in_tuple, 0);
+  
+} // BaseCase_()
+
+
+
+
+void npt::SingleBandwidthAlg::DepthFirstRecursion_(std::vector<SingleNode*>& nodes) {
+  
+  // check for symmetry and pruning
+  bool can_prune = CheckNodeList_(nodes);
+  
+  if (can_prune) {
+    
+    // note that this will count prunes based on symmetry too
+    num_prunes_++;
+    
+  }
+  else {
+    
+    // look over all the nodes, see if they are leaves, and if not, which one 
+    // to splits
+    bool all_leaves = nodes[0]->is_leaf();
+    index_t split_index = 0;
+    // if node 0 is not a leaf, then use it, otherwise don't
+    index_t split_count = all_leaves ? -1 : nodes[0]->count();
+
+    // loop over the other nodes, check for leaves and who to split
+    for (index_t i = 1; i < tuple_size_; i++) {
+     
+      if (!(nodes[i]->is_leaf())) {
+        
+        all_leaves = false;
+        
+        if (nodes[i]->count() > split_count) {
+          split_count = nodes[i]->count();
+          split_index = i;
+        }
+        
+      } // not a leaf
+      
+    } // for i
+    
+    // Do we do the base case or recurse?
+    
+    if (all_leaves) {
+      
+      BaseCase_(nodes);
+      
+    } 
+    else {
+      
+      SingleNode* split_node = nodes[split_index];
+      
+      // left child first
+      nodes[split_index] = split_node->left();
+      DepthFirstRecursion_(nodes);
+      
+      // now do the other child
+      nodes[split_index] = split_node->right();
+      DepthFirstRecursion_(nodes);
+      
+      // This is important for calls above this one.
+      nodes[split_index] = split_node;
+      
+      
+    } // not a base case
+          
+  } // can't prune
+  
+} // DepthFirstRecursion_
@@ -29,6 +29,7 @@ namespace npt {
    // input params
    index_t num_points_;
    index_t tuple_size_;
+    index_t leaf_size_;
    
    // the matcher
    SingleMatcher matcher_;
@@ -41,10 +42,78 @@ namespace npt {
    // the number of times we pruned a tuple
    int num_prunes_;
    
-    NPointNode* tree_;
+    // define the tree
+    typedef BinarySpaceTree<DHrectBound<2>, arma::mat> SingleNode; 
+    
+    arma::Col<index_t> old_from_new_index_;
+    
+    SingleNode* tree_;
+    
+    
+    ////////////////////// functions /////////////////////////
+    
+    bool CheckNodeList_(std::vector<SingleNode*>& nodes);
+    
+    void BaseCaseHelper_(std::vector<std::vector<index_t> >& point_sets,
+                         std::vector<bool>& permutation_ok,
+                         std::vector<index_t>& points_in_tuple,
+                         int k);
+    
+    void BaseCase_(std::vector<SingleNode*>& nodes);
+    
+    void DepthFirstRecursion_(std::vector<SingleNode*>& nodes);
+    
    
  public:
    
+    /**
+     * Requires the matcher bounds, data, parameters
+     */
+    SingleBandwidthAlg(arma::mat& data, arma::colvec weights, index_t n,
+                       index_t leaf_size,
+                       arma::mat& lower_bds, arma::mat& upper_bds) : 
+                       matcher_(n, lower_bds, upper_bds)
+    {
+      
+      data_points_ = data;
+      data_weights_ = weights;
+      num_points_ = data_points_.n_cols;
+      tuple_size_ = n;
+      leaf_size_ = leaf_size;
+      
+      num_tuples_ = 0;
+      weighted_num_tuples_ = 0.0;
+      num_prunes_ = 0;
+      
+      tree_ = tree::MakeKdTreeMidpoint<SingleNode, double>(data_points_, 
+                                                           leaf_size_, 
+                                                           old_from_new_index_);
+        
+      // IMPORTANT: need to permute the weights here
+      
+      
+      
+    } // constructor
+    
+    int num_tuples() {
+      return num_tuples_;
+    } 
+    
+    double weighted_num_tuples() {
+      return weighted_num_tuples_;
+    }
+    
+    
+    /**
+     * Actually run the algorithm.
+     */
+    void ComputeCounts() {
+            
+      std::vector<SingleNode*> node_list(tuple_size_, tree_);
+      
+      DepthFirstRecursion_(node_list);
+      
+    } // ComputeCounts()
    
    
    
@@ -75,6 +75,10 @@ namespace npt {
                       index_t tuple_ind_1, index_t tuple_ind_2,
                       std::vector<bool>& permutation_ok);
    
+    index_t num_permutations() {
+      return num_permutations_;
+    }
+    
    
  }; // class