source: trunk/src/vector.cc @ 132

Last change on this file since 132 was 132, checked in by Peter, 18 years ago

modified stream constructor to ignore empty lines

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.2 KB
Line 
1// $Id: vector.cc 132 2004-08-09 13:24:38Z peter $
2
3#include "vector.h"
4#include "random_singleton.h"
5
6#include <iostream>
7#include <sstream>
8#include <string>
9#include <vector>
10#include <utility>
11
12
13
14
15
16namespace theplu {
17namespace gslapi {
18
19
20
21  vector::vector(void)
22    : v_(NULL)
23  {
24  }
25
26
27
28  vector::vector(const size_t n, const double init_value)
29  {
30    v_ = gsl_vector_alloc(n);
31    set_all(init_value);
32  }
33
34
35
36  vector::vector(const vector& other, const std::vector<size_t>& index)
37  {
38    if (index.size()){
39      v_ = gsl_vector_alloc(index.size());
40      for (size_t i=0; i<index.size(); i++) 
41        gsl_vector_set(v_, i, other(index[i]));
42    }
43    else
44      v_ = other.gsl_vector_copy();
45  }
46
47
48
49  vector::vector(gsl_vector* vec)
50    : v_(vec)
51  {
52  }
53
54
55
56  vector::vector(std::istream& is)
57  {
58    using namespace std;
59 
60    // read the data file and store in stl vectors (dynamically expandable)
61    std::vector<std::vector<double> > data_matrix;
62    u_int nof_columns=0;
63    u_int nof_rows = 0;
64    string s;
65    for (nof_rows = 0; getline(is, s, '\n'); nof_rows++) {
66      istringstream line(s);
67      std::vector<double> v;
68      string tmp_string;
69      while (line >> tmp_string) {
70        if(!is.good()) {
71          // Jari. change to throw exception, remove unnecessary includes
72          cerr << "vector::vector(std::istream&): "
73               << "error reading data file!" << endl;
74          exit(1);
75        }
76        double t=atof(tmp_string.c_str());
77        // Here we should check that it actually was a double!!!!!
78        v.push_back(t);
79      }
80
81      // Ignoring empty lines
82      if(!v.size()){
83        nof_rows--;
84        continue;
85      }
86
87      if(nof_columns==0)
88        nof_columns=v.size();
89      else if(v.size()!=nof_columns) {
90        // Jari. change to throw exception, remove unnecessary includes
91        cerr << "vector::vector(std::istream&) data file error: "
92             << "line" << nof_rows+1 << " has " << v.size() 
93             << " columns; expected " << nof_columns
94             << " columns"
95             << endl; 
96        exit(1);
97      }
98      data_matrix.push_back(v);
99    }
100 
101    // manipulate the state of the stream to be good
102    is.clear(std::ios::goodbit);
103 
104    // convert the data to a gsl vector and check that data file is a
105    // column vector or a row vector
106    if(nof_columns==1) {
107      v_ = gsl_vector_alloc ( nof_rows );
108      for(u_int i=0;i<nof_rows;i++) 
109        gsl_vector_set( v_, i, data_matrix[i][0] );
110    }
111    else if(nof_rows==1){
112      v_ = gsl_vector_alloc ( nof_columns );
113      for(u_int i=0;i<nof_columns;i++) 
114        gsl_vector_set( v_, i, data_matrix[0][i] );
115    }
116    else {
117      // Jari. change to throw exception, remove unnecessary includes
118      cerr << "vector::vector(std::istream&) data file error: "
119           << "file has " << nof_rows << " rows and " << nof_columns
120           << " columns; expected a row vector or a column vector"
121           << endl;
122      exit(1);
123    }
124  }
125
126
127
128  vector::~vector(void)
129  {
130    if (v_) {
131      gsl_vector_free(v_);
132      v_=NULL;
133    }
134  }
135
136
137
138  gsl_vector* vector::gsl_vector_copy(void) const
139  {
140    gsl_vector* vec = gsl_vector_alloc(size());
141    gsl_vector_memcpy(vec,v_);
142    return vec;
143  }
144
145
146
147
148
149  std::pair<size_t,size_t> vector::minmax_index(void) const
150  {
151    size_t min_index=0;
152    size_t max_index=0;
153    void gsl_vector_minmax_index (const gsl_vector * v_, size_t * 
154                                  min_index, size_t * max_index);
155    return std::pair<size_t,size_t>(min_index, max_index);
156  }
157
158
159
160
161  std::pair<size_t,size_t> vector::minmax_index(const std::vector<size_t>& subset ) const
162  {
163    size_t min_index = subset[0];
164    size_t max_index = subset[0];
165    for (unsigned int i=1; i<subset.size(); i++) {
166      if (gsl_vector_get( v_, subset[i]) < gsl_vector_get( v_, min_index))
167        min_index = subset[i];
168      else if (gsl_vector_get( v_, subset[i]) > gsl_vector_get( v_, max_index))
169        max_index = subset[i];
170    }
171    return std::pair<size_t,size_t>(min_index, max_index);
172  }
173
174
175
176  vector vector::mul_elements( const vector& other ) const
177  {
178    vector res( *this );
179    gsl_vector_mul( res.v_, other.v_ );
180    return res;
181  }
182
183  void vector::shuffle(void) const
184  {
185    vector vec(*this);
186    theplu::cpptools::random_singleton* 
187      rnd=theplu::cpptools::random_singleton::get_instance();
188    for (size_t i=0; i<vec.size(); i++){
189      size_t index = rnd->get_uniform_int(vec.size()-i);
190      gsl_vector_set(v_,i,vec(index));
191      vec(index)=vec(vec.size()-i-1);
192    }
193  }
194
195  double vector::sum(void) const
196  {
197    double sum = 0;
198    for (u_int i=0; i<size(); i++)
199      sum += gsl_vector_get( v_, i );
200    return( sum );
201  } 
202
203
204
205  vector vector::operator-(void) const
206  {
207    vector res( *this );
208    gsl_vector_scale (res.v_,-1.);
209    return res;
210  }
211
212
213
214  double vector::operator*( const vector &other ) const
215  {
216    // Jari, check for gsl_support
217    double res = 0.0;;
218    for ( size_t i = 0; i < size(); ++i ) 
219      res += other[i] * (*this)[i];
220    return res;
221  }
222
223
224
225  vector vector::operator+(const vector &other) const
226  {
227    vector res(*this);
228    gsl_vector_add(res.v_,other.v_);
229    return res;
230  }
231
232
233
234  vector vector::operator-( const vector &other ) const
235  {
236    vector res( *this );
237    gsl_vector_sub(res.v_,other.v_);
238    return res;
239  }
240
241
242
243  vector& vector::operator=( const vector& other )
244  {
245    if( this != &other ) {
246      if ( v_ )
247        gsl_vector_free( v_ );
248      v_ = other.gsl_vector_copy();
249    }
250    return *this;
251  } 
252
253
254
255  std::ostream& theplu::gslapi::operator<<(std::ostream& s, const vector& a)
256  {
257    s.setf(std::ios::fixed);
258
259    for (size_t j = 0; j < a.size(); ++j) {
260      s << a[j];
261      if ( (j+1)<a.size() )
262        s << " ";
263    }
264
265    return s;
266  }
267
268
269}} // of namespace gslapi and namespace thep
Note: See TracBrowser for help on using the repository browser.