source: trunk/yat/classifier/MatrixLookupWeighted.cc @ 1486

Last change on this file since 1486 was 1486, checked in by Jari Häkkinen, 13 years ago

Addresses #436.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date ID
File size: 9.2 KB
Line 
1// $Id$
2
3/*
4  Copyright (C) 2006 Jari Häkkinen, Peter Johansson, Markus Ringnér
5  Copyright (C) 2007 Jari Häkkinen, Peter Johansson
6  Copyright (C) 2008 Peter Johansson
7
8  This file is part of the yat library, http://dev.thep.lu.se/yat
9
10  The yat library is free software; you can redistribute it and/or
11  modify it under the terms of the GNU General Public License as
12  published by the Free Software Foundation; either version 3 of the
13  License, or (at your option) any later version.
14
15  The yat library is distributed in the hope that it will be useful,
16  but WITHOUT ANY WARRANTY; without even the implied warranty of
17  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18  General Public License for more details.
19
20  You should have received a copy of the GNU General Public License
21  along with this program; if not, write to the Free Software
22  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
23  02111-1307, USA.
24*/
25
26#include "MatrixLookupWeighted.h"
27#include "MatrixLookup.h"
28#include "yat/utility/DataIterator.h"
29#include "yat/utility/Matrix.h"
30#include "yat/utility/MatrixWeighted.h"
31#include "yat/utility/WeightIterator.h"
32
33#include <algorithm>
34#include <cassert>
35#include <fstream>
36#include <vector>
37
38namespace theplu {
39namespace yat {
40namespace classifier {
41
42  MatrixLookupWeighted::MatrixLookupWeighted(const utility::MatrixWeighted& m,
43                                             const utility::Index& rows,
44                                             const utility::Index& columns)
45    : column_index_(columns), row_index_(rows)
46  {
47    utility::Matrix* data = new utility::Matrix(m.rows(), m.columns());
48    utility::Matrix* weight = new utility::Matrix(m.rows(), m.columns());
49    for (size_t i=0; i<m.rows(); ++i)
50      for (size_t j=0; j<m.columns(); ++j) {
51        (*data)(i,j) = m(i,j).data();
52        (*weight)(i,j) = m(i,j).weight();
53      }
54    // smart pointers are taking ownership
55    data_ = MatrixP(data);
56    weights_ = MatrixP(weight);
57  }
58
59
60  MatrixLookupWeighted::MatrixLookupWeighted(const utility::MatrixWeighted& m)
61    : column_index_(utility::Index(m.columns())), 
62      row_index_(utility::Index(m.rows()))
63  {
64    utility::Matrix* data = new utility::Matrix(m.rows(), m.columns());
65    utility::Matrix* weight = new utility::Matrix(m.rows(), m.columns());
66    for (size_t i=0; i<m.rows(); ++i)
67      for (size_t j=0; j<m.columns(); ++j) {
68        (*data)(i,j) = m(i,j).data();
69        (*weight)(i,j) = m(i,j).weight();
70      } 
71  // smart pointers are taking ownership
72    data_ = MatrixP(data);
73    weights_ = MatrixP(weight);
74  }
75
76
77  MatrixLookupWeighted::MatrixLookupWeighted(const utility::Matrix& data, 
78                                             const utility::Matrix& weights,
79                                             const bool own)
80    : data_(MatrixP(&data, own)), weights_(MatrixP(&weights, own)) 
81  {
82    assert(data.rows()==weights.rows());
83    assert(data.columns()==weights.columns());
84    row_index_ = utility::Index(data.rows());
85    column_index_ = utility::Index(data.columns());
86    assert(validate());
87  }
88
89
90  MatrixLookupWeighted::MatrixLookupWeighted(const MatrixLookup& ml)
91    : column_index_(ml.column_index_), data_(ml.data_),
92      row_index_(ml.row_index_), 
93      weights_(MatrixP(new utility::Matrix(data_->rows(),data_->columns(),1.0)))
94  { 
95    assert(validate());
96  }
97 
98
99  MatrixLookupWeighted::MatrixLookupWeighted(const utility::Matrix& data, 
100                                             const utility::Matrix& weights, 
101                                             const utility::Index& row, 
102                                             const utility::Index& col)
103    : column_index_(col), data_(MatrixP(&data, false)), 
104      row_index_(row), weights_(MatrixP(&weights, false))
105  {
106    assert(validate());
107  }
108 
109
110
111  /*
112  MatrixLookupWeighted::MatrixLookupWeighted(const utility::Matrix& data,
113                                             const utility::Matrix& weights,
114                                             const utility::Index& index,
115                                             const bool row)
116    : data_(MatrixP(new utility::Matrix(data), false)),
117      weights_(MatrixP(new utility::Matrix(weights), false))
118  {
119    assert(data.rows()==weights.rows());
120    assert(data.columns()==weights.columns());
121    if (row){
122      row_index_=index;
123      column_index_ = utility::Index(data.columns());
124    }
125    else{
126      column_index_=index;
127      row_index_ = utility::Index(data.rows());
128    }
129    assert(validate());
130  }
131  */
132
133  /*
134  MatrixLookupWeighted::MatrixLookupWeighted(const MatrixLookup& dv,
135                                             const MatrixLookup& wv)
136    : DataLookup2D(dv), data_(dv.data_), weights_(dv.data_)
137  {
138  }
139  */
140
141
142  MatrixLookupWeighted::MatrixLookupWeighted(const MatrixLookupWeighted& other)
143    : column_index_(other.column_index_), data_(other.data_), 
144      row_index_(other.row_index_), weights_(other.weights_)
145  {
146    assert(validate());
147  }
148
149
150
151  MatrixLookupWeighted::MatrixLookupWeighted(const MatrixLookupWeighted& other,
152                                             const utility::Index& row, 
153                                             const utility::Index& col)
154    : column_index_(utility::Index(other.column_index_, col)), 
155      data_(other.data_), row_index_(utility::Index(other.row_index_, row)), 
156      weights_(other.weights_)
157  {
158    assert(validate());
159  }
160 
161
162
163  MatrixLookupWeighted::MatrixLookupWeighted(const MatrixLookupWeighted& other, 
164                                             const utility::Index& index, 
165                                             bool row)
166    : data_(other.data_), 
167      weights_(other.weights_)
168  {
169    if (row){
170      row_index_ = utility::Index(other.row_index_, index);
171      column_index_= other.column_index_;
172    }
173    else{
174      column_index_ = utility::Index(other.column_index_, index);
175      row_index_= other.row_index_;
176    }
177    assert(validate());
178  }
179 
180
181
182  MatrixLookupWeighted::MatrixLookupWeighted(const size_t rows, 
183                                             const size_t columns, 
184                                             const double value,
185                                             const double weight)
186    : column_index_(utility::Index(std::vector<size_t>(rows, 0))),
187      data_(MatrixP(new utility::Matrix(1,1,value))),
188      row_index_(utility::Index(std::vector<size_t>(columns, 0))), 
189      weights_(MatrixP(new utility::Matrix(1,1,weight)))
190  {
191    assert(validate());
192  }
193
194 
195  MatrixLookupWeighted::MatrixLookupWeighted(std::istream& is, char sep)
196    : data_(MatrixP(new utility::Matrix(is,sep)))
197  {
198    column_index_ = utility::Index(data_->columns());
199    row_index_ = utility::Index(data_->rows());
200    utility::Matrix weights;
201    utility::nan(*data_,weights);
202    // Peter, should be possible to avoid this copying
203    weights_= MatrixP(new utility::Matrix(weights));
204    assert(validate());
205  }
206 
207
208  MatrixLookupWeighted::~MatrixLookupWeighted(void)
209  {
210  }
211
212
213
214  MatrixLookupWeighted::const_iterator MatrixLookupWeighted::begin(void) const
215  {
216    return const_iterator(const_iterator::iterator_type(*this, 0, 0), 1);
217  }
218
219
220  MatrixLookupWeighted::const_column_iterator
221  MatrixLookupWeighted::begin_column(size_t i) const
222  {
223    return const_column_iterator(const_column_iterator::iterator_type(*this,0, 
224                                                                      i),
225                                 columns());
226  }
227
228
229  MatrixLookupWeighted::const_iterator
230  MatrixLookupWeighted::begin_row(size_t i) const
231  {
232    return const_row_iterator(const_row_iterator::iterator_type(*this,i, 0), 1);
233  }
234
235
236  size_t MatrixLookupWeighted::columns(void) const
237  {
238    return column_index_.size();
239  }
240
241
242  double MatrixLookupWeighted::data(size_t row, size_t column) const
243  {
244    assert(row<rows());
245    assert(column<columns());
246    assert(row_index_[row]<data_->rows());
247    assert(column_index_[column]<data_->columns());
248    return (*data_)(row_index_[row], column_index_[column]);
249  }
250
251
252
253  MatrixLookupWeighted::const_iterator MatrixLookupWeighted::end(void) const
254  {
255    return const_iterator(const_iterator::iterator_type(*this, rows(), 0), 1);
256  }
257
258
259  MatrixLookupWeighted::const_iterator
260  MatrixLookupWeighted::end_column(size_t i) const
261  {
262    return const_column_iterator(const_column_iterator::iterator_type(*this, 
263                                                                      rows(),i),
264                                 columns());
265  }
266
267
268  MatrixLookupWeighted::const_row_iterator
269  MatrixLookupWeighted::end_row(size_t i) const
270  {
271    return const_row_iterator(const_row_iterator::iterator_type(*this,i+1,0),1);
272  }
273
274
275  size_t MatrixLookupWeighted::rows(void) const
276  {
277    return row_index_.size();
278  }
279
280
281  bool MatrixLookupWeighted::validate(void) const
282  {
283    for (size_t i=0; i<row_index_.size(); ++i)
284      if (row_index_[i]>=data_->rows()) {
285        std::cerr << i << " " << row_index_[i] << std::endl;
286        return false;
287      }
288    for (size_t i=0; i<column_index_.size(); ++i)
289      if (column_index_[i]>=data_->columns())
290        return false;
291    for (size_t i=0; i<row_index_.size(); ++i)
292      if (row_index_[i]>=weights_->rows())
293        return false;
294    for (size_t i=0; i<column_index_.size(); ++i)
295      if (column_index_[i]>=weights_->columns())
296        return false;
297    return true;
298  }
299
300
301  double MatrixLookupWeighted::weight(size_t row, size_t column) const
302  {
303    assert(row<rows());
304    assert(column<columns());
305    assert(row_index_[row]<weights_->rows());
306    assert(column_index_[column]<weights_->columns());
307    return (*weights_)(row_index_[row], column_index_[column]);
308  }
309
310
311
312  bool MatrixLookupWeighted::weighted(void) const 
313  {
314    return true;
315  }
316
317
318
319  std::pair<double, double> 
320  MatrixLookupWeighted::operator()(const size_t row, const size_t column) const
321  { 
322    return std::make_pair(data(row, column), weight(row,column));
323  }
324
325
326
327  const MatrixLookupWeighted& MatrixLookupWeighted::operator=
328  (const MatrixLookupWeighted& other)
329  {
330    if (this!=&other){
331      column_index_=other.column_index_;
332      row_index_=other.row_index_;
333      data_ = other.data_;
334      weights_ = other.weights_;
335    }
336    assert(validate());
337    return *this;
338  }
339
340
341}}} // of namespace classifier, yat, and theplu
Note: See TracBrowser for help on using the repository browser.