LCOV - code coverage report
Current view: top level - routines - multidimensional_scaling.hpp (source / functions) Hit Total Coverage
Test: clean.info Lines: 53 53 100.0 %
Date: 2013-05-24 Functions: 9 32 28.1 %
Branches: 63 323 19.5 %

           Branch data     Line data    Source code
       1                 :            : /* This software is distributed under BSD 3-clause license (see LICENSE file).
       2                 :            :  *
       3                 :            :  * Copyright (c) 2012-2013 Sergey Lisitsyn
       4                 :            :  */
       5                 :            : 
       6                 :            : #ifndef TAPKEE_MDS_H_
       7                 :            : #define TAPKEE_MDS_H_
       8                 :            : 
       9                 :            : /* Tapkee includes */
      10                 :            : #include <tapkee/defines.hpp>
      11                 :            : #include <tapkee/utils/time.hpp>
      12                 :            : /* End of Tapkee includes */
      13                 :            : 
      14                 :            : namespace tapkee
      15                 :            : {
      16                 :            : namespace tapkee_internal
      17                 :            : {
      18                 :            : 
      19                 :            : template <class RandomAccessIterator>
      20                 :          2 : Landmarks select_landmarks_random(RandomAccessIterator begin, RandomAccessIterator end, ScalarType ratio)
      21                 :            : {
      22                 :          2 :         Landmarks landmarks;
      23   [ +  -  #  # ]:          2 :         landmarks.reserve(end-begin);
      24 [ +  + ][ #  # ]:        102 :         for (RandomAccessIterator iter=begin; iter!=end; ++iter)
      25 [ +  - ][ #  # ]:        100 :                 landmarks.push_back(iter-begin);
      26 [ +  - ][ +  - ]:          2 :         tapkee::random_shuffle(landmarks.begin(),landmarks.end());
         [ +  - ][ #  # ]
         [ #  # ][ #  # ]
      27 [ +  - ][ +  - ]:          2 :         landmarks.erase(landmarks.begin() + static_cast<IndexType>(landmarks.size()*ratio),landmarks.end());
         [ +  - ][ #  # ]
         [ #  # ][ #  # ]
      28                 :          2 :         return landmarks;
      29                 :            : }
      30                 :            : 
      31                 :            : template <class RandomAccessIterator, class PairwiseCallback>
      32                 :          1 : DenseSymmetricMatrix compute_distance_matrix(RandomAccessIterator begin, RandomAccessIterator /*end*/,
      33                 :            :                                              const Landmarks& landmarks, PairwiseCallback callback)
      34                 :            : {
      35 [ +  - ][ +  - ]:          1 :         timed_context context("Multidimensional scaling distance matrix computation");
         [ +  - ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
      36                 :            : 
      37                 :          1 :         const IndexType n_landmarks = landmarks.size();
      38   [ +  -  #  #  :          1 :         DenseSymmetricMatrix distance_matrix(n_landmarks,n_landmarks);
             #  #  #  # ]
      39                 :            : 
      40                 :            : #pragma omp parallel shared(begin,landmarks,distance_matrix,callback) default(none)
      41                 :            :         {
      42                 :            :                 IndexType i_index_iter,j_index_iter;
      43 [ +  + ][ #  # ]:          5 : #pragma omp for nowait
         [ #  # ][ #  # ]
      44                 :            :                 for (i_index_iter=0; i_index_iter<n_landmarks; ++i_index_iter)
      45                 :            :                 {
      46 [ +  + ][ +  + ]:        341 :                         for (j_index_iter=i_index_iter; j_index_iter<n_landmarks; ++j_index_iter)
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
      47                 :            :                         {
      48                 :        314 :                                 ScalarType d = callback.distance(begin[landmarks[i_index_iter]],begin[landmarks[j_index_iter]]);
      49                 :        287 :                                 d *= d;
      50                 :        287 :                                 distance_matrix(i_index_iter,j_index_iter) = d;
      51                 :        313 :                                 distance_matrix(j_index_iter,i_index_iter) = d;
      52                 :            :                         }
      53                 :            :                 }
      54                 :            :         }
      55                 :          1 :         return distance_matrix;
      56                 :            : }
      57                 :            : 
      58                 :            : template <class RandomAccessIterator, class PairwiseCallback>
      59                 :          1 : DenseMatrix triangulate(RandomAccessIterator begin, RandomAccessIterator end, PairwiseCallback distance_callback,
      60                 :            :                         const Landmarks& landmarks, const DenseVector& landmark_distances_squared, 
      61                 :            :                         EigendecompositionResult& landmarks_embedding, IndexType target_dimension)
      62                 :            : {
      63 [ +  - ][ +  - ]:          1 :         timed_context context("Landmark triangulation");
         [ +  - ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
      64                 :            :         
      65                 :          1 :         const IndexType n_vectors = end-begin;
      66                 :          1 :         const IndexType n_landmarks = landmarks.size();
      67                 :            : 
      68         [ +  - ]:         51 :         bool* to_process = new bool[n_vectors];
           [ +  +  #  # ]
           [ #  #  #  # ]
           [ #  #  #  # ]
                 [ #  # ]
      69 [ +  - ][ #  # ]:          1 :         std::fill(to_process,to_process+n_vectors,true);
         [ #  # ][ #  # ]
      70                 :            :         
      71 [ +  - ][ #  # ]:          1 :         DenseMatrix embedding(n_vectors,target_dimension);
         [ #  # ][ #  # ]
      72                 :            : 
      73 [ +  + ][ #  # ]:         26 :         for (IndexType index_iter=0; index_iter<n_landmarks; ++index_iter)
         [ #  # ][ #  # ]
      74                 :            :         {
      75                 :         25 :                 to_process[landmarks[index_iter]] = false;
      76 [ +  - ][ +  - ]:         25 :                 embedding.row(landmarks[index_iter]).noalias() = landmarks_embedding.first.row(index_iter);
                 [ +  - ]
           [ +  -  #  # ]
         [ #  # ][ #  # ]
           [ #  #  #  # ]
         [ #  # ][ #  # ]
           [ #  #  #  # ]
         [ #  # ][ #  # ]
                 [ #  # ]
      77                 :            :         }
      78                 :            : 
      79 [ +  + ][ #  # ]:          3 :         for (IndexType i=0; i<target_dimension; ++i)
         [ #  # ][ #  # ]
      80 [ +  - ][ +  - ]:          2 :                 landmarks_embedding.first.col(i).array() /= landmarks_embedding.second(i);
         [ +  - ][ +  - ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
      81                 :            : 
      82                 :            : #pragma omp parallel shared(begin,end,to_process,distance_callback,landmarks, \
      83                 :            :                 landmarks_embedding,landmark_distances_squared,embedding) default(none)
      84                 :            :         {
      85                 :          3 :                 DenseVector distances_to_landmarks(n_landmarks);
      86                 :            :                 IndexType index_iter;
      87 [ +  - ][ #  # ]:          3 : #pragma omp for nowait
         [ #  # ][ #  # ]
      88                 :            :                 for (index_iter=0; index_iter<n_vectors; ++index_iter)
      89                 :            :                 {
      90 [ +  + ][ +  + ]:         55 :                         if (!to_process[index_iter])
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
      91                 :         27 :                                 continue;
      92                 :            : 
      93 [ +  + ][ #  # ]:        516 :                         for (IndexType i=0; i<n_landmarks; ++i)
         [ #  # ][ #  # ]
      94                 :            :                         {
      95                 :        497 :                                 ScalarType d = distance_callback.distance(begin[index_iter],begin[landmarks[i]]);
      96                 :        505 :                                 distances_to_landmarks(i) = d*d;
      97                 :            :                         }
      98                 :            :                         //distances_to_landmarks.array().square();
      99                 :            : 
     100                 :         19 :                         distances_to_landmarks -= landmark_distances_squared;
     101                 :         30 :                         embedding.row(index_iter).noalias() = -0.5*landmarks_embedding.first.transpose()*distances_to_landmarks;
     102                 :            :                 }
     103                 :            :         }
     104                 :            : 
     105 [ +  - ][ #  # ]:          1 :         delete[] to_process;
         [ #  # ][ #  # ]
     106                 :            : 
     107                 :          1 :         return embedding;
     108                 :            : }
     109                 :            : 
     110                 :            : template <class RandomAccessIterator, class PairwiseCallback>
     111                 :         10 : DenseSymmetricMatrix compute_distance_matrix(RandomAccessIterator begin, RandomAccessIterator end, 
     112                 :            :                                              PairwiseCallback callback)
     113                 :            : {
     114 [ +  - ][ +  - ]:         10 :         timed_context context("Multidimensional scaling distance matrix computation");
         [ +  - ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
     115                 :            : 
     116                 :         10 :         const IndexType n_vectors = end-begin;
     117   [ +  +  #  #  :         10 :         DenseSymmetricMatrix distance_matrix(n_vectors,n_vectors);
             #  #  #  # ]
     118                 :            : 
     119                 :            : #pragma omp parallel shared(begin,distance_matrix,callback) default(none)
     120                 :            :         {
     121                 :            :                 IndexType i_index_iter,j_index_iter;
     122 [ +  # ][ #  # ]:         68 : #pragma omp for nowait
         [ #  # ][ +  + ]
     123                 :            :                 for (i_index_iter=0; i_index_iter<n_vectors; ++i_index_iter)
     124                 :            :                 {
     125 [ +  + ][ +  + ]:       7335 :                         for (j_index_iter=i_index_iter; j_index_iter<n_vectors; ++j_index_iter)
         [ #  # ][ #  # ]
         [ #  # ][ #  # ]
            [ +  + ][ + ]
     126                 :            :                         {
     127                 :       7019 :                                 ScalarType d = callback.distance(begin[i_index_iter],begin[j_index_iter]);
     128                 :       6927 :                                 d *= d;
     129                 :       6927 :                                 distance_matrix(i_index_iter,j_index_iter) = d;
     130                 :       7032 :                                 distance_matrix(j_index_iter,i_index_iter) = d;
     131                 :            :                         }
     132                 :            :                 }
     133                 :            :         }
     134                 :         10 :         return distance_matrix;
     135                 :            : }
     136                 :            : 
     137                 :            : } // End of namespace tapkee_internal
     138                 :            : } // End of namespace tapkee
     139                 :            : 
     140                 :            : #endif

Generated by: LCOV version 1.9