Browse Source

bbw tutorial example

Former-commit-id: 66a3d3a8fa270fbf9f453b13f1926144670397b7
Alec Jacobson 11 years ago
parent
commit
df835027b4

+ 49 - 0
include/igl/deform_skeleton.cpp

@@ -0,0 +1,49 @@
+#include "deform_skeleton.h"
+void igl::deform_skeleton(
+  const Eigen::MatrixXd & C,
+  const Eigen::MatrixXi & BE,
+  const std::vector<
+    Eigen::Affine3d,Eigen::aligned_allocator<Eigen::Affine3d> > & vA,
+  Eigen::MatrixXd & CT,
+  Eigen::MatrixXi & BET)
+{
+  using namespace Eigen;
+  assert(BE.rows() == (int)vA.size());
+  CT.resize(2*BE.rows(),C.cols());
+  BET.resize(BE.rows(),2);
+  for(int e = 0;e<BE.rows();e++)
+  {
+    BET(e,0) = 2*e; 
+    BET(e,1) = 2*e+1;
+    Affine3d a = vA[e];
+    Vector3d c0 = C.row(BE(e,0));
+    Vector3d c1 = C.row(BE(e,1));
+    CT.row(2*e) =   a * c0;
+    CT.row(2*e+1) = a * c1;
+  }
+
+}
+
+IGL_INLINE void igl::deform_skeleton(
+  const Eigen::MatrixXd & C,
+  const Eigen::MatrixXi & BE,
+  const Eigen::MatrixXd & T,
+  Eigen::MatrixXd & CT,
+  Eigen::MatrixXi & BET)
+{
+  using namespace Eigen;
+  assert(BE.rows() == (int)vA.size());
+  CT.resize(2*BE.rows(),C.cols());
+  BET.resize(BE.rows(),2);
+  for(int e = 0;e<BE.rows();e++)
+  {
+    BET(e,0) = 2*e; 
+    BET(e,1) = 2*e+1;
+    Affine3d a;
+    a.matrix() = T.block(e*4,0,4,3).transpose();
+    Vector3d c0 = C.row(BE(e,0));
+    Vector3d c1 = C.row(BE(e,1));
+    CT.row(2*e) =   a * c0;
+    CT.row(2*e+1) = a * c1;
+  }
+}

+ 39 - 0
include/igl/deform_skeleton.h

@@ -0,0 +1,39 @@
+#ifndef IGL_DEFORM_SKELETON_H
+#define IGL_DEFORM_SKELETON_H
+#include "igl_inline.h"
+#include <Eigen/Core>
+#include <Eigen/Geometry>
+#include <vector>
+namespace igl
+{
+  // Deform a skeleton.
+  //
+  // Inputs:
+  //   C  #C by 3 list of joint positions
+  //   BE  #BE by 2 list of bone edge indices
+  //   vA  #BE list of bone transformations
+  // Outputs
+  //   CT  #BE*2 by 3 list of deformed joint positions
+  //   BET  #BE by 2 list of bone edge indices (maintains order)
+  //
+  IGL_INLINE void deform_skeleton(
+    const Eigen::MatrixXd & C,
+    const Eigen::MatrixXi & BE,
+    const std::vector<
+      Eigen::Affine3d,Eigen::aligned_allocator<Eigen::Affine3d> > & vA,
+    Eigen::MatrixXd & CT,
+    Eigen::MatrixXi & BET);
+  // Inputs:
+  //   T  #BE*4 by 3 list of stacked transformation matrix
+  IGL_INLINE void deform_skeleton(
+    const Eigen::MatrixXd & C,
+    const Eigen::MatrixXi & BE,
+    const Eigen::MatrixXd & T,
+    Eigen::MatrixXd & CT,
+    Eigen::MatrixXi & BET);
+}
+  
+#ifndef IGL_STATIC_LIBRARY
+#  include "deform_skeleton.cpp"
+#endif
+#endif

+ 2 - 1
include/igl/viewer/TODOs.txt

@@ -1,5 +1,5 @@
+- data.lines, data.points should not concatenate colors with coordinates
 - snap to canonical recenters origin but trackball does not
-+ snap to canonical view key shortcut is not working
 - rewrite in libigl style
 - separate various class into their own .h/.cpp pairs
 - remove use of double underscores (http://stackoverflow.com/a/224420/148668)
@@ -14,6 +14,7 @@
 - zoom with pan rather than scaling
 - refresh draw while resizing
 - use constructor initializer list rather than complicated constructor
++ snap to canonical view key shortcut is not working
 + resize TwBar with window
 + trackball should be able to drag over TwBar
 + don't zoom on horizontal scale

+ 1 - 1
include/igl/viewer/Viewer.cpp.REMOVED.git-id

@@ -1 +1 @@
-5c850c8fd21c8bb05e0ccf7fb0dbcc8e9eeeb6b8
+b02c6fbf6cb464fdaf3bfb95194978587d5cbd50

+ 8 - 0
include/igl/viewer/Viewer.h

@@ -371,6 +371,14 @@ namespace igl
                       const Eigen::Matrix<char,Eigen::Dynamic,Eigen::Dynamic>& B);
 
     void add_points(const Eigen::MatrixXd& P,  const Eigen::MatrixXd& C);
+    // Sets edges given a list of edge vertices and edge indices. In constrast
+    // to `add_edges` this will (purposefully) clober existing edges.
+    //
+    // Inputs:
+    //   P  #P by 3 list of vertex positions
+    //   E  #E by 2 list of edge indices into P
+    //   C  #E|1 by 3 color(s)
+    void set_edges (const Eigen::MatrixXd& P, const Eigen::MatrixXi& E, const Eigen::MatrixXd& C);
     void add_edges (const Eigen::MatrixXd& P1, const Eigen::MatrixXd& P2, const Eigen::MatrixXd& C);
     void add_label (const Eigen::VectorXd& P,  const std::string& str);
 

+ 0 - 3
tutorial/401_BiharmonicDeformation/main.cpp

@@ -1,8 +1,5 @@
 #include <igl/colon.h>
-#include <igl/cotmatrix.h>
 #include <igl/harmonic.h>
-#include <igl/jet.h>
-#include <igl/min_quad_with_fixed.h>
 #include <igl/readOBJ.h>
 #include <igl/readDMAT.h>
 #include <igl/viewer/Viewer.h>

+ 0 - 3
tutorial/402_PolyharmonicDeformation/main.cpp

@@ -1,8 +1,5 @@
 #include <igl/colon.h>
-#include <igl/cotmatrix.h>
 #include <igl/harmonic.h>
-#include <igl/jet.h>
-#include <igl/min_quad_with_fixed.h>
 #include <igl/readOBJ.h>
 #include <igl/viewer/Viewer.h>
 #include <algorithm>

+ 11 - 0
tutorial/403_BoundedBiharmonicWeights/CMakeLists.txt

@@ -0,0 +1,11 @@
+cmake_minimum_required(VERSION 2.6)
+project(403_BoundedBiharmonicWeights)
+
+include("../CMakeLists.shared")
+
+set(SOURCES
+${PROJECT_SOURCE_DIR}/main.cpp
+)
+
+add_executable(${PROJECT_NAME}_bin ${SOURCES} ${SHARED_SOURCES})
+target_link_libraries(${PROJECT_NAME}_bin ${SHARED_LIBRARIES})

+ 169 - 0
tutorial/403_BoundedBiharmonicWeights/main.cpp

@@ -0,0 +1,169 @@
+// Don't use static library for this example because of Mosek complications
+//#define IGL_NO_MOSEK
+#ifdef IGL_NO_MOSEK
+#undef IGL_STATIC_LIBRARY
+#endif
+#include <igl/boundary_conditions.h>
+#include <igl/colon.h>
+#include <igl/column_to_quats.h>
+#include <igl/forward_kinematics.h>
+#include <igl/jet.h>
+#include <igl/lbs_matrix.h>
+#include <igl/deform_skeleton.h>
+#include <igl/normalize_row_sums.h>
+#include <igl/readDMAT.h>
+#include <igl/readMESH.h>
+#include <igl/readTGF.h>
+#include <igl/viewer/Viewer.h>
+#include <igl/bbw/bbw.h>
+
+#include <Eigen/Geometry>
+#include <Eigen/StdVector>
+#include <vector>
+#include <algorithm>
+#include <iostream>
+
+typedef 
+  std::vector<Eigen::Quaterniond,Eigen::aligned_allocator<Eigen::Quaterniond> >
+  RotationList;
+
+const Eigen::RowVector3d sea_green(70./255.,252./255.,167./255.);
+int selected = 0;
+Eigen::MatrixXd V,W,U,C,M;
+Eigen::MatrixXi T,F,BE;
+Eigen::VectorXi P;
+RotationList pose;
+double anim_t = 1.0;
+double anim_t_dir = -0.03;
+
+bool pre_draw(igl::Viewer & viewer)
+{
+  using namespace Eigen;
+  using namespace std;
+  if(viewer.options.is_animating)
+  {
+    // Interpolate pose and identity
+    RotationList anim_pose(pose.size());
+    for(int e = 0;e<pose.size();e++)
+    {
+      anim_pose[e] = pose[e].slerp(anim_t,Quaterniond::Identity());
+    }
+    // Propogate relative rotations via FK to retrieve absolute transformations
+    RotationList vQ;
+    vector<Vector3d> vT;
+    igl::forward_kinematics(C,BE,P,anim_pose,vQ,vT);
+    const int dim = C.cols();
+    MatrixXd T(BE.rows()*(dim+1),dim);
+    for(int e = 0;e<BE.rows();e++)
+    {
+      Affine3d a = Affine3d::Identity();
+      a.translate(vT[e]);
+      a.rotate(vQ[e]);
+      T.block(e*(dim+1),0,dim+1,dim) =
+        a.matrix().transpose().block(0,0,dim+1,dim);
+    }
+    // Compute deformation via LBS as matrix multiplication
+    U = M*T;
+
+    // Also deform skeleton edges
+    MatrixXd CT;
+    MatrixXi BET;
+    igl::deform_skeleton(C,BE,T,CT,BET);
+    
+    viewer.set_vertices(U);
+    viewer.set_edges(CT,BET,sea_green);
+    viewer.compute_normals();
+    anim_t += anim_t_dir;
+    anim_t_dir *= (anim_t>=1.0 || anim_t<=0.0?-1.0:1.0);
+  }
+  return false;
+}
+
+void set_color(igl::Viewer &viewer)
+{
+  Eigen::MatrixXd C;
+  igl::jet(W.col(selected).eval(),true,C);
+  viewer.set_colors(C);
+}
+
+bool key_down(igl::Viewer &viewer, unsigned char key, int mods)
+{
+  switch(key)
+  {
+    case ' ':
+      viewer.options.is_animating = !viewer.options.is_animating;
+      break;
+    case '.':
+      selected++;
+      selected = std::min(std::max(selected,0),(int)W.cols()-1);
+      set_color(viewer);
+      break;
+    case ',':
+      selected--;
+      selected = std::min(std::max(selected,0),(int)W.cols()-1);
+      set_color(viewer);
+      break;
+  }
+}
+
+int main(int argc, char *argv[])
+{
+  using namespace Eigen;
+  using namespace std;
+  igl::readMESH("../shared/hand.mesh",V,T,F);
+  U=V;
+  igl::readTGF("../shared/hand.tgf",C,BE);
+  // retrieve parents for forward kinematics
+  P.resize(BE.rows(),1);
+  for(int e = 0;e<BE.rows();e++)
+  {
+    P(e) = -1;
+    for(int f = 0;f<BE.rows();f++)
+    {
+      if(BE(e,0) == BE(f,1))
+      {
+        P(e) = f;
+      }
+    }
+  }
+
+  // Read pose as matrix of quaternions per row
+  MatrixXd Q;
+  igl::readDMAT("../shared/hand-pose.dmat",Q);
+  igl::column_to_quats(Q,pose);
+  assert(pose.size() == BE.rows());
+
+  // List of boundary indices (aka fixed value indices into VV)
+  VectorXi b;
+  // List of boundary conditions of each weight function
+  MatrixXd bc;
+  igl::boundary_conditions(V,T,C,VectorXi(),BE,MatrixXi(),b,bc);
+
+  // compute BBW weights matrix
+  igl::BBWData bbw_data;
+  // only a few iterations for sake of demo
+  bbw_data.active_set_params.max_iter = 8;
+  if(!igl::bbw(V,T,b,bc,bbw_data,W))
+  {
+    return false;
+  }
+  // Normalize weights to sum to one
+  igl::normalize_row_sums(W,W);
+  // precompute linear blend skinning matrix
+  igl::lbs_matrix(V,W,M);
+
+  // Plot the mesh with pseudocolors
+  igl::Viewer viewer;
+  viewer.set_mesh(U, F);
+  set_color(viewer);
+  viewer.set_edges(C,BE,sea_green);
+  viewer.options.show_lines = false;
+  viewer.options.show_overlay_depth = false;
+  viewer.options.line_width = 1;
+  viewer.options.trackball_angle.normalize();
+  viewer.callback_pre_draw = &pre_draw;
+  viewer.callback_key_down = &key_down;
+  viewer.options.is_animating = false;
+  viewer.options.animation_max_fps = 30.;
+  viewer.launch();
+}

+ 18 - 0
tutorial/cmake/FindLIBIGL.cmake

@@ -42,6 +42,24 @@ if(LIBIGL_USE_STATIC_LIBRARY)
   endif(NOT LIBIGL_LIBRARY)
   set(LIBIGL_LIBRARIES ${LIBIGL_LIBRARIES}  ${LIBIGL_LIBRARY})
 
+  FIND_LIBRARY( LIBIGLBBW_LIBRARY NAMES iglbbw PATHS ${LIBIGL_LIB_DIRS})
+  if(NOT LIBIGLBBW_LIBRARY)
+    set(LIBIGL_FOUND FALSE)
+  endif(NOT LIBIGLBBW_LIBRARY)
+  set(LIBIGL_LIBRARIES ${LIBIGL_LIBRARIES}  ${LIBIGLBBW_LIBRARY})
+  FIND_LIBRARY( LIBIGLMOSEK_LIBRARY NAMES iglmosek PATHS ${LIBIGL_LIB_DIRS})
+  if(NOT LIBIGLMOSEK_LIBRARY)
+    set(LIBIGL_FOUND FALSE)
+  endif(NOT LIBIGLMOSEK_LIBRARY)
+  set(LIBIGL_LIBRARIES ${LIBIGL_LIBRARIES}  ${LIBIGLMOSEK_LIBRARY})
+  find_package(Mosek REQUIRED)
+  if(MOSEK_FOUND)
+    set(LIBIGL_INCLUDE_DIR ${LIBIGL_INCLUDE_DIR}  ${MOSEK_INCLUDE_DIR})
+    set(LIBIGL_LIBRARIES ${LIBIGL_LIBRARIES}  ${MOSEK_LIBRARIES})
+  else(MOSEK_FOUND)
+    set(LIBIGL_FOUND FALSE)
+  endif(MOSEK_FOUND)
+
   FIND_LIBRARY( LIBIGLEMBREE_LIBRARY NAMES iglembree PATHS ${LIBIGL_LIB_DIRS})
   if(NOT LIBIGLEMBREE_LIBRARY)
     set(LIBIGL_FOUND FALSE)

+ 1 - 0
tutorial/images/hand-bbw.jpg.REMOVED.git-id

@@ -0,0 +1 @@
+65841362faf22ff5e76c1375fc3ee6d62f2a78fe

+ 65 - 0
tutorial/readme.md

@@ -862,11 +862,76 @@ igl::harmonic(V,F,b,bc,k,Z);
 ![The `PolyharmonicDeformation` example deforms a flat domain (left) into a bump as a
 solution to various $k$-harmonic PDEs.](images/bump-k-harmonic.jpg)
 
+## Bounded Biharmonic Weights
+In computer animation, shape deformation is often referred to as "skinning".
+Constraints are posed as relative rotations of internal rigid "bones" inside a
+character. The deformation method, or skinning method, determines how the
+surface of the character (i.e. its skin) should move as a function of the bone
+rotations.
+
+The most popular technique is linear blend skinning. Each point on the shape
+computes its new location as a linear combination of bone transformations:
+
+ $\mathbf{x}' = \sum\limits_{i = 0}^m w_i(\mathbf{x}) \mathbf{T}_i
+ \left(\begin{array}{c}\mathbf{x}_i\\1\end{array}\right),$
+
+where $w_i(\mathbf{x})$ is the scalar _weight function_ of the ith bone evaluated at
+$\mathbf{x}$ and $\mathbf{T}_i$ is the bone transformation as a $4 \times 3$
+matrix.
+
+This formula is embarassingly parallel (computation at one point does not
+depend on shared data need by computation at another point). It is often
+implemented as a vertex shader. The weights and rest positions for each vertex
+are sent as vertex shader _attribtues_ and bone transformations are sent as
+_uniforms_. Then vertices are transformed within the vertex shader, just in
+time for rendering.
+
+As the skinning formula is linear (hence its name), we can write it as matrix
+multiplication:
+
+ $\mathbf{X}' = \mathbf{M} \mathbf{T},$
+
+where $\mathbx{X}'$ is $n \times 3$ stack of deformed positions as row
+vectors, $\mathbf{M}$ is a $n \times m\cdot dim$ matrix containing weights and
+rest positions and $\mathbf{T}$ is a $m\cdot (dim+1) \times dim$ stack of
+transposed bone transformations.
+
+Traditionally, the weight functions $w_j$ are painted manually by skilled
+rigging professionals. Modern techniques now exist to compute weight functions
+automatically given the shape and a description of the skeleton (or in general
+any handle structure such as a cage, collection of points, selected regions,
+etc.).
+
+Bounded biharmonic weights are one such technique that casts weight computation
+as a constrained optimization problem [#jacobson_2011][]. The weights enforce
+smoothness by minimizing a smoothness energy: the familiar Laplacian energy:
+
+ $\sum\limits_{i = 0}^m \int_S (\Delta w_i)^2 dA$
+  
+subject to constraints which enforce interpolation of handle constraints:
+
+ $w_i(\mathbf{x}) = \begin{cases} 1 & \text{ if } \mathbf{x} \in H_i\\ 0 & \text{ otherwise }
+ \end{cases},$
+
+where $H_i$ is the ith handle, and constraints which enforce non-negativity,
+parition of unity and encourage sparsity:
+
+ $0\le w_i \le 1$ and $\sum\limits_{i=0}^m w_i = 1.$
+
+This is a quadratic programming problem and libigl solves it using its active
+set solver or by calling out to Mosek.
+
+![The example `BoundedBiharmonicWeights` computes weights for a tetrahedral
+mesh given a skeleton (top) and then animates a linear blend skinning
+deformation (bottom).](images/hand-bbw.jpg)
+
 [#botsch_2004]: Matrio Botsch and Leif Kobbelt. "An Intuitive Framework for
 Real-Time Freeform Modeling," 2004.
 [#jacobson_thesis_2013]: Alec Jacobson,
 _Algorithms and Interfaces for Real-Time Deformation of 2D and 3D Shapes_,
 2013.
+[#jacobson_2011]: Alec Jacobson, Ilya Baran, Jovan Popović, and Olga Sorkin.
+"Bounded Biharmonic Weights for Real-Time Deformation," 2011.
 [#jacobson_mixed_2010]: Alec Jacobson, Elif Tosun, Olga Sorkine, and Denis
 Zorin. "Mixed Finite Elements for Variational Surface Modeling," 2010.
 [#kazhdan_2012]: Michael Kazhdan, Jake Solomon, Mirela Ben-Chen,

+ 81 - 0
tutorial/shared/hand-pose.dmat

@@ -0,0 +1,81 @@
+1 80
+0.61618142029453848
+0.21865691494778211
+-0.074398815808914262
+0.75297704283434619
+0.016137740405538198
+0.0075736234381567954
+-0.01978957096844574
+0.9996452302909804
+-0.0031783575088926955
+-0.0021852197161076415
+0.025195604356894566
+0.99967509940951804
+0
+0
+0
+1
+0
+0
+0
+1
+-0.06502916547759717
+-0.018228928610557554
+-0.10960238233656835
+0.99167849204525316
+0.65862930444592327
+-0.28749500451272225
+0.047433167985789754
+0.69376087831482691
+0.71088878530987565
+0.13051734182532121
+-0.024300471136119696
+0.69066044153843431
+0.12918716086559792
+0.42797711290327217
+-0.23378675140515992
+0.86341764121767584
+0.5777562073061252
+-0.10739726497891923
+0.027027719666595052
+0.80866129792634012
+0
+0
+0
+1
+0.25064723245034182
+0.062362300611793375
+0.055537180373327777
+0.96447007725671363
+0.68233863926855265
+-0.043115978196253735
+0.066162231411292738
+0.72675824929614585
+0.53643420171330392
+-0.34234055129349178
+0.079151746448423155
+0.76731759735246274
+0.37119350486218949
+-0.20651629016779105
+0.015792865442698228
+0.90516130564935249
+0.46007394537647717
+-0.033329312116362406
+-0.011156066109515889
+0.88718468422788521
+0.36118819983774919
+0.046197896403048229
+0.048746682979448855
+0.93007128735619049
+-0.080434792811141881
+0.24444948358806326
+-0.37172646406082138
+0.89196083433960183
+0.34016846232137932
+0.054859354927957418
+0.046210032868129636
+0.93762492569317746
+0.45250460737446119
+0.27010308530876809
+0.1405545502553627
+0.83816962603883571

+ 1 - 0
tutorial/shared/hand.mesh.REMOVED.git-id

@@ -0,0 +1 @@
+cbe696a1280075895d31178317b30d722e0af9f3

+ 43 - 0
tutorial/shared/hand.tgf

@@ -0,0 +1,43 @@
+   1 0.11429762 -0.76590296 -0.11140733          0          1          0          1          1          1          0          0          0          1          1  1.0061325 0.88859267          0    1 
+   2 0.10546931 -0.57479244 -0.084241846 2.7982368e+199 3.7150555e-27 6.0174446e+175          1          1          1          0          0          0          1 0.99850723   1.182572 0.91575815          0    2 
+   3 -0.24048184 -0.34974589 0.083588056 2.0750757e-322 3.9525252e-322 2.5296161e-321          1          1          1          0          0          0          1 0.65901122  1.2530116  1.0835881          0    4 
+   4 -0.41528797 -0.19724228 0.15108035 6.953262e-310 6.953262e-310 6.953262e-310          1          1          1          0          0          0          1 0.86263359  1.1364105  1.1510803          0    6 
+   5 -0.53251027 -0.043704366   0.220108          0          1          0          1          1          1          0          0          0          1 0.91920809  1.1497897   1.220108          0    8 
+   6 -0.16766723 0.12476223 -0.049004442 1.0869444e-322 3.9525252e-322 2.5296161e-321          1          1          1          0          0          0          1 0.72621948  1.7247165 0.95099556          0   10 
+   7 -0.22466126 0.29596172 -0.045066203 6.9532642e-310 6.9532642e-310 6.9532642e-310          1          1          1          0          0          0          1 0.92992065  1.1895161  0.9549338          0   12 
+   8 -0.29358606 0.48060868 -0.018623317 5.4347221e-323 5.4347221e-323 3.9525252e-323          1          1          1          0          0          0          1 0.93107519  1.1835472 0.98137668          0   14 
+   9 -0.3288664 0.60375622 0.0023204123 5.4347221e-323 5.4347221e-323 3.9525252e-323          1          1          1          0          0          0          1 0.94905516  1.1441107  1.0023204          0   16 
+  10 -0.0045596564 0.17429162 -0.063580889 5.4347221e-323 5.4347221e-323 3.9525252e-323          1          1          1          0          0          0          1 0.87811436  1.7574269 0.93641911          0   18 
+  11 -0.035737249 0.36523184 -0.049670244          0          1          0          1          1          1          0          0          0          1 0.95065421  1.2251482 0.95032976          0   20 
+  12 -0.063819738 0.56050302 -0.026654403 5.4347221e-323 5.4347221e-323 3.9525252e-323          1          1          1          0          0          0          1 0.95901298  1.2354098  0.9733456          0   22 
+  13 -0.080813321 0.71934417 0.0051323801 6.9532642e-310 6.9532642e-310 6.9532642e-310          1          1          1          0          0          0          1 0.96724796  1.2232573  1.0051324          0   24 
+  14 0.15872682  0.1411217 -0.041958815          0          1          0          1          1          1          0          0          0          1   1.044204  1.7326665 0.95804118          0   26 
+  15 0.16999343 0.35280711 -0.0202035          0          1          0          1          1          1          0          0          0          1  1.0246357  1.2351488  0.9797965          0   28 
+  16  0.1877683 0.51268791 0.010202099 6.9533558e-309 6.9532426e-310          0          1          1          1          0          0          0          1   1.017182  1.1791852  1.0102021          0   30 
+  17 0.20129211 0.64288912 0.029292092 6.9398703e-310 6.9398703e-310 6.9398703e-310          1          1          1          0          0          0          1   1.002674  1.1159674  1.0292921          0   32 
+  18 0.34208039 0.052370865          0 6.9398715e-310 2.3912777e-321          0          1          1          1          0          0          0          1  1.2359671  1.6383093          1          0   34 
+  19 0.40904104 0.18734965 0.027992567          0          1          0          1          1          1          0          0          0          1   1.081477  1.1621961  1.0279926          0   36 
+  20 0.49154459 0.32275416 0.058527871 6.9398703e-310 6.9398703e-310 6.9398703e-310          1          1          1          0          0          0          1   1.074279  1.1322133  1.0585279          0   38 
+  21 0.54469073 0.41146288 0.082016739 6.9398703e-310 2.3912777e-321          0          1          1          1          0          0          0          1  1.0377236  1.0807633  1.0820167          0   40 
+#
+   6    7
+   2    6
+   2   18
+   2   10
+   2   14
+   2    3
+  18   19
+  10   11
+   4    5
+  16   17
+   1    2
+   8    9
+  14   15
+  20   21
+  19   20
+  15   16
+  11   12
+   3    4
+  12   13
+   7    8
+#