00001
00002 #ifndef vnl_matrix_fixed_txx_
00003 #define vnl_matrix_fixed_txx_
00004
00005
00006 #include "vnl_matrix_fixed.h"
00007
00008 #include <vcl_cmath.h>
00009 #include <vcl_iostream.h>
00010 #include <vcl_cstdlib.h>
00011 #include <vcl_cassert.h>
00012
00013 #include <vnl/vnl_error.h>
00014 #include <vnl/vnl_math.h>
00015 #include <vnl/vnl_vector_fixed.h>
00016
00017 template<class T, unsigned nrows, unsigned ncols>
00018 void
00019 vnl_matrix_fixed<T,nrows,ncols>::add( const T* a, const T* b, T* r )
00020 {
00021 unsigned int count = nrows*ncols;
00022 while ( count-- )
00023 *(r++) = *(a++) + *(b++);
00024 }
00025
00026
00027 template<class T, unsigned nrows, unsigned ncols>
00028 void
00029 vnl_matrix_fixed<T,nrows,ncols>::add( const T* a, T b, T* r )
00030 {
00031 unsigned int count = nrows*ncols;
00032 while ( count-- )
00033 *(r++) = *(a++) + b;
00034 }
00035
00036 template<class T, unsigned nrows, unsigned ncols>
00037 void
00038 vnl_matrix_fixed<T,nrows,ncols>::sub( const T* a, const T* b, T* r )
00039 {
00040 unsigned int count = nrows*ncols;
00041 while ( count-- )
00042 *(r++) = *(a++) - *(b++);
00043 }
00044
00045 template<class T, unsigned nrows, unsigned ncols>
00046 void
00047 vnl_matrix_fixed<T,nrows,ncols>::sub( const T* a, T b, T* r )
00048 {
00049 unsigned int count = nrows*ncols;
00050 while ( count-- )
00051 *(r++) = *(a++) - b;
00052 }
00053
00054 template<class T, unsigned nrows, unsigned ncols>
00055 void
00056 vnl_matrix_fixed<T,nrows,ncols>::sub( T a, const T* b, T* r )
00057 {
00058 unsigned int count = nrows*ncols;
00059 while ( count-- )
00060 *(r++) = a - *(b++);
00061 }
00062
00063 template<class T, unsigned nrows, unsigned ncols>
00064 void
00065 vnl_matrix_fixed<T,nrows,ncols>::mul( const T* a, const T* b, T* r )
00066 {
00067 unsigned int count = nrows*ncols;
00068 while ( count-- )
00069 *(r++) = *(a++) * *(b++);
00070 }
00071
00072 template<class T, unsigned nrows, unsigned ncols>
00073 void
00074 vnl_matrix_fixed<T,nrows,ncols>::mul( const T* a, T b, T* r )
00075 {
00076 unsigned int count = nrows*ncols;
00077 while ( count-- )
00078 *(r++) = *(a++) * b;
00079 }
00080
00081 template<class T, unsigned nrows, unsigned ncols>
00082 void
00083 vnl_matrix_fixed<T,nrows,ncols>::div( const T* a, const T* b, T* r )
00084 {
00085 unsigned int count = nrows*ncols;
00086 while ( count-- )
00087 *(r++) = *(a++) / *(b++);
00088 }
00089
00090 template<class T, unsigned nrows, unsigned ncols>
00091 void
00092 vnl_matrix_fixed<T,nrows,ncols>::div( const T* a, T b, T* r )
00093 {
00094 unsigned int count = nrows*ncols;
00095 while ( count-- )
00096 *(r++) = *(a++) / b;
00097 }
00098
00099 template<class T, unsigned nrows, unsigned ncols>
00100 bool
00101 vnl_matrix_fixed<T,nrows,ncols>::equal( const T* a, const T* b )
00102 {
00103 unsigned int count = nrows*ncols;
00104 while ( count-- )
00105 if ( *(a++) != *(b++) ) return false;
00106 return true;
00107 }
00108
00109
00110
00111
00112 template<class T, unsigned nrows, unsigned ncols>
00113 vnl_matrix_fixed<T,nrows,ncols>&
00114 vnl_matrix_fixed<T,nrows,ncols>::fill (T value)
00115 {
00116 for (unsigned int i = 0; i < nrows; ++i)
00117 for (unsigned int j = 0; j < ncols; ++j)
00118 this->data_[i][j] = value;
00119 return *this;
00120 }
00121
00122
00123 template<class T, unsigned nrows, unsigned ncols>
00124 vnl_matrix_fixed<T,nrows,ncols>&
00125 vnl_matrix_fixed<T,nrows,ncols>::fill_diagonal (T value)
00126 {
00127 for (unsigned int i = 0; i < nrows && i < ncols; ++i)
00128 this->data_[i][i] = value;
00129 return *this;
00130 }
00131
00132
00133 template<class T, unsigned nrows, unsigned ncols>
00134 vnl_matrix_fixed<T,nrows,ncols>&
00135 vnl_matrix_fixed<T,nrows,ncols>::set_diagonal(vnl_vector<T> const& diag)
00136 {
00137 assert(diag.size() >= nrows || diag.size() >= ncols);
00138
00139
00140
00141 for (unsigned int i = 0; i < nrows && i < ncols; ++i)
00142 this->data_[i][i] = diag[i];
00143 return *this;
00144 }
00145
00146
00147 template<class T, unsigned nrows, unsigned ncols>
00148 void
00149 vnl_matrix_fixed<T,nrows,ncols>::print(vcl_ostream& os) const
00150 {
00151 for (unsigned int i = 0; i < nrows; ++i)
00152 {
00153 os << this->data_[i][0];
00154 for (unsigned int j = 1; j < ncols; ++j)
00155 os << ' ' << this->data_[i][j];
00156 os << '\n';
00157 }
00158 }
00159
00160
00161 template <class T, unsigned nrows, unsigned ncols>
00162 vnl_matrix_fixed<T,nrows,ncols>
00163 vnl_matrix_fixed<T,nrows,ncols>::apply(T (*f)(T const&)) const
00164 {
00165 vnl_matrix_fixed<T,nrows,ncols> ret;
00166 vnl_c_vector<T>::apply(this->data_[0], rows()*cols(), f, ret.data_block());
00167 return ret;
00168 }
00169
00170 template <class T, unsigned nrows, unsigned ncols>
00171 vnl_matrix_fixed<T,nrows,ncols>
00172 vnl_matrix_fixed<T,nrows,ncols>::apply(T (*f)(T)) const
00173 {
00174 vnl_matrix_fixed<T,nrows,ncols> ret;
00175 vnl_c_vector<T>::apply(this->data_[0], rows()*cols(), f, ret.data_block());
00176 return ret;
00177 }
00178
00179
00180
00181
00182 template<class T, unsigned nrows, unsigned ncols>
00183 vnl_matrix_fixed<T,ncols,nrows>
00184 vnl_matrix_fixed<T,nrows,ncols>::transpose() const
00185 {
00186 vnl_matrix_fixed<T,ncols,nrows> result;
00187 for (unsigned int i = 0; i < cols(); ++i)
00188 for (unsigned int j = 0; j < rows(); ++j)
00189 result(i,j) = this->data_[j][i];
00190 return result;
00191 }
00192
00193 template<class T, unsigned nrows, unsigned ncols>
00194 vnl_matrix_fixed<T,ncols,nrows>
00195 vnl_matrix_fixed<T,nrows,ncols>::conjugate_transpose() const
00196 {
00197 vnl_matrix_fixed<T,ncols,nrows> result(transpose());
00198 vnl_c_vector<T>::conjugate(result.begin(),
00199 result.begin(),
00200 result.size());
00201 return result;
00202 }
00203
00204 template<class T, unsigned nrows, unsigned ncols>
00205 vnl_matrix_fixed<T,nrows,ncols>&
00206 vnl_matrix_fixed<T,nrows,ncols>::update (vnl_matrix<T> const& m,
00207 unsigned top, unsigned left)
00208 {
00209 const unsigned int bottom = top + m.rows();
00210 const unsigned int right = left + m.cols();
00211 #ifndef NDEBUG
00212 if (nrows < bottom || ncols < right)
00213 vnl_error_matrix_dimension ("update",
00214 bottom, right, m.rows(), m.cols());
00215 #endif
00216 for (unsigned int i = top; i < bottom; ++i)
00217 for (unsigned int j = left; j < right; ++j)
00218 this->data_[i][j] = m(i-top,j-left);
00219 return *this;
00220 }
00221
00222
00223 template<class T, unsigned nrows, unsigned ncols>
00224 vnl_matrix<T>
00225 vnl_matrix_fixed<T,nrows,ncols>::extract (unsigned rowz, unsigned colz,
00226 unsigned top, unsigned left) const
00227 {
00228 vnl_matrix<T> result(rowz, colz);
00229 this->extract( result, top, left );
00230 return result;
00231 }
00232
00233
00234 template<class T, unsigned nrows, unsigned ncols>
00235 void
00236 vnl_matrix_fixed<T,nrows,ncols>::extract (vnl_matrix<T>& sub_matrix,
00237 unsigned top, unsigned left) const
00238 {
00239 unsigned int rowz = sub_matrix.rows();
00240 unsigned int colz = sub_matrix.cols();
00241 #ifndef NDEBUG
00242 unsigned int bottom = top + rowz;
00243 unsigned int right = left + colz;
00244 if ((nrows < bottom) || (ncols < right))
00245 vnl_error_matrix_dimension ("extract",
00246 nrows, ncols, bottom, right);
00247 #endif
00248 for (unsigned int i = 0; i < rowz; ++i)
00249 for (unsigned int j = 0; j < colz; ++j)
00250 sub_matrix(i,j) = this->data_[top+i][left+j];
00251 }
00252
00253
00254 template<class T, unsigned nrows, unsigned ncols>
00255 vnl_matrix_fixed<T,nrows,ncols>&
00256 vnl_matrix_fixed<T,nrows,ncols>::copy_in(T const *p)
00257 {
00258 T* dp = this->data_block();
00259 unsigned int i = nrows * ncols;
00260 while (i--)
00261 *dp++ = *p++;
00262 return *this;
00263 }
00264
00265 template<class T, unsigned nrows, unsigned ncols>
00266 void vnl_matrix_fixed<T,nrows,ncols>::copy_out(T *p) const
00267 {
00268 T const* dp = this->data_block();
00269 unsigned int i = nrows*ncols;
00270 while (i--)
00271 *p++ = *dp++;
00272 }
00273
00274 template<class T, unsigned nrows, unsigned ncols>
00275 vnl_matrix_fixed<T,nrows,ncols>&
00276 vnl_matrix_fixed<T,nrows,ncols>::set_identity()
00277 {
00278
00279
00280 for (unsigned int i = 0; i < nrows; ++i)
00281 for (unsigned int j = 0; j < ncols; ++j)
00282 this->data_[i][j] = T(0);
00283 for (unsigned int i = 0; i < nrows && i < ncols; ++i)
00284 this->data_[i][i] = T(1);
00285 return *this;
00286 }
00287
00288
00289
00290 template<class T, unsigned nrows, unsigned ncols>
00291 vnl_matrix_fixed<T,nrows,ncols>&
00292 vnl_matrix_fixed<T,nrows,ncols>::normalize_rows()
00293 {
00294 for (unsigned int i = 0; i < nrows; ++i)
00295 {
00296 abs_t norm(0);
00297 for (unsigned int j = 0; j < ncols; ++j)
00298 norm += vnl_math_squared_magnitude( this->data_[i][j] );
00299
00300 if (norm != 0)
00301 {
00302 typedef typename vnl_numeric_traits<abs_t>::real_t real_t;
00303 real_t scale = real_t(1)/vcl_sqrt((real_t)norm);
00304 for (unsigned int j = 0; j < ncols; ++j)
00305 {
00306
00307
00308 this->data_[i][j] *= T(scale);
00309 }
00310 }
00311 }
00312 return *this;
00313 }
00314
00315 template<class T, unsigned nrows, unsigned ncols>
00316 vnl_matrix_fixed<T,nrows,ncols>&
00317 vnl_matrix_fixed<T,nrows,ncols>::normalize_columns()
00318 {
00319 for (unsigned int j = 0; j < ncols; ++j) {
00320 abs_t norm(0);
00321 for (unsigned int i = 0; i < nrows; ++i)
00322 norm += vnl_math_squared_magnitude( this->data_[i][j] );
00323
00324 if (norm != 0)
00325 {
00326 typedef typename vnl_numeric_traits<abs_t>::real_t real_t;
00327 real_t scale = real_t(1)/vcl_sqrt((real_t)norm);
00328 for (unsigned int i = 0; i < nrows; ++i)
00329 {
00330
00331
00332 this->data_[i][j] *= T(scale);
00333 }
00334 }
00335 }
00336 return *this;
00337 }
00338
00339 template<class T, unsigned nrows, unsigned ncols>
00340 vnl_matrix_fixed<T,nrows,ncols>&
00341 vnl_matrix_fixed<T,nrows,ncols>::scale_row(unsigned row_index, T value)
00342 {
00343 #ifndef NDEBUG
00344 if (row_index >= nrows)
00345 vnl_error_matrix_row_index("scale_row", row_index);
00346 #endif
00347 for (unsigned int j = 0; j < ncols; ++j)
00348 this->data_[row_index][j] *= value;
00349 return *this;
00350 }
00351
00352 template<class T, unsigned nrows, unsigned ncols>
00353 vnl_matrix_fixed<T,nrows,ncols>&
00354 vnl_matrix_fixed<T,nrows,ncols>::scale_column(unsigned column_index, T value)
00355 {
00356 #ifndef NDEBUG
00357 if (column_index >= ncols)
00358 vnl_error_matrix_col_index("scale_column", column_index);
00359 #endif
00360 for (unsigned int j = 0; j < nrows; ++j)
00361 this->data_[j][column_index] *= value;
00362 return *this;
00363 }
00364
00365
00366 template<class T, unsigned nrows, unsigned ncols>
00367 vnl_matrix<T>
00368 vnl_matrix_fixed<T,nrows,ncols>::get_n_rows (unsigned row, unsigned n) const
00369 {
00370 #ifndef NDEBUG
00371 if (row + n > nrows)
00372 vnl_error_matrix_row_index ("get_n_rows", row);
00373 #endif
00374
00375
00376 return vnl_matrix<T>(data_[row], n, ncols);
00377 }
00378
00379 template<class T, unsigned nrows, unsigned ncols>
00380 vnl_matrix<T>
00381 vnl_matrix_fixed<T,nrows,ncols>::get_n_columns (unsigned column, unsigned n) const
00382 {
00383 #ifndef NDEBUG
00384 if (column + n > ncols)
00385 vnl_error_matrix_col_index ("get_n_columns", column);
00386 #endif
00387
00388 vnl_matrix<T> result(nrows, n);
00389 for (unsigned int c = 0; c < n; ++c)
00390 for (unsigned int r = 0; r < nrows; ++r)
00391 result(r, c) = this->data_[r][column + c];
00392 return result;
00393 }
00394
00395
00396 template<class T, unsigned nrows, unsigned ncols>
00397 vnl_vector_fixed<T,ncols> vnl_matrix_fixed<T,nrows,ncols>::get_row(unsigned row_index) const
00398 {
00399 #ifdef ERROR_CHECKING
00400 if (row_index >= nrows)
00401 vnl_error_matrix_row_index ("get_row", row_index);
00402 #endif
00403
00404 vnl_vector_fixed<T,ncols> v;
00405 for (unsigned int j = 0; j < ncols; ++j)
00406 v[j] = this->data_[row_index][j];
00407 return v;
00408 }
00409
00410
00411 template<class T, unsigned nrows, unsigned ncols>
00412 vnl_vector_fixed<T,nrows> vnl_matrix_fixed<T,nrows,ncols>::get_column(unsigned column_index) const
00413 {
00414 #ifdef ERROR_CHECKING
00415 if (column_index >= ncols)
00416 vnl_error_matrix_col_index ("get_column", column_index);
00417 #endif
00418
00419 vnl_vector_fixed<T,nrows> v;
00420 for (unsigned int j = 0; j < nrows; ++j)
00421 v[j] = this->data_[j][column_index];
00422 return v;
00423 }
00424
00425
00426 template<class T, unsigned nrows, unsigned ncols>
00427 vnl_vector<T> vnl_matrix_fixed<T,nrows,ncols>::get_diagonal() const
00428 {
00429 vnl_vector<T> v(nrows < ncols ? nrows : ncols);
00430 for (unsigned int j = 0; j < nrows && j < ncols; ++j)
00431 v[j] = this->data_[j][j];
00432 return v;
00433 }
00434
00435
00436
00437 template<class T, unsigned nrows, unsigned ncols>
00438 vnl_matrix_fixed<T,nrows,ncols>&
00439 vnl_matrix_fixed<T,nrows,ncols>::set_row(unsigned row_index, T const *v)
00440 {
00441 for (unsigned int j = 0; j < ncols; ++j)
00442 this->data_[row_index][j] = v[j];
00443 return *this;
00444 }
00445
00446 template<class T, unsigned nrows, unsigned ncols>
00447 vnl_matrix_fixed<T,nrows,ncols>&
00448 vnl_matrix_fixed<T,nrows,ncols>::set_row(unsigned row_index, vnl_vector<T> const &v)
00449 {
00450 if (v.size() >= ncols)
00451 set_row(row_index,v.data_block());
00452 else
00453 for (unsigned int j = 0; j < v.size(); ++j)
00454 this->data_[row_index][j] = v[j];
00455 return *this;
00456 }
00457
00458 template<class T, unsigned nrows, unsigned ncols>
00459 vnl_matrix_fixed<T,nrows,ncols>&
00460 vnl_matrix_fixed<T,nrows,ncols>::set_row(unsigned row_index, vnl_vector_fixed<T,ncols> const &v)
00461 {
00462 set_row(row_index,v.data_block());
00463 return *this;
00464 }
00465
00466 template<class T, unsigned nrows, unsigned ncols>
00467 vnl_matrix_fixed<T,nrows,ncols>&
00468 vnl_matrix_fixed<T,nrows,ncols>::set_row(unsigned row_index, T v)
00469 {
00470 for (unsigned int j = 0; j < ncols; ++j)
00471 this->data_[row_index][j] = v;
00472 return *this;
00473 }
00474
00475
00476
00477 template<class T, unsigned nrows, unsigned ncols>
00478 vnl_matrix_fixed<T,nrows,ncols>&
00479 vnl_matrix_fixed<T,nrows,ncols>::set_column(unsigned column_index, T const *v)
00480 {
00481 for (unsigned int i = 0; i < nrows; ++i)
00482 this->data_[i][column_index] = v[i];
00483 return *this;
00484 }
00485
00486 template<class T, unsigned nrows, unsigned ncols>
00487 vnl_matrix_fixed<T,nrows,ncols>&
00488 vnl_matrix_fixed<T,nrows,ncols>::set_column(unsigned column_index, vnl_vector<T> const &v)
00489 {
00490 if (v.size() >= nrows)
00491 set_column(column_index,v.data_block());
00492 else
00493 for (unsigned int i = 0; i < v.size(); ++i)
00494 this->data_[i][column_index] = v[i];
00495 return *this;
00496 }
00497
00498 template<class T, unsigned nrows, unsigned ncols>
00499 vnl_matrix_fixed<T,nrows,ncols>&
00500 vnl_matrix_fixed<T,nrows,ncols>::set_column(unsigned column_index, vnl_vector_fixed<T,nrows> const &v)
00501 {
00502 set_column(column_index,v.data_block());
00503 return *this;
00504 }
00505
00506 template<class T, unsigned nrows, unsigned ncols>
00507 vnl_matrix_fixed<T,nrows,ncols>&
00508 vnl_matrix_fixed<T,nrows,ncols>::set_column(unsigned column_index, T v)
00509 {
00510 for (unsigned int j = 0; j < nrows; ++j)
00511 this->data_[j][column_index] = v;
00512 return *this;
00513 }
00514
00515
00516 template<class T, unsigned nrows, unsigned ncols>
00517 vnl_matrix_fixed<T,nrows,ncols>&
00518 vnl_matrix_fixed<T,nrows,ncols>::set_columns(unsigned starting_column, vnl_matrix<T> const& m)
00519 {
00520 for (unsigned int j = 0; j < m.cols() && starting_column+j < ncols; ++j)
00521 for (unsigned int i = 0; i < nrows && i < m.rows(); ++i)
00522 this->data_[i][starting_column + j] = m(i,j);
00523 return *this;
00524 }
00525
00526
00527 template <class T, unsigned nrows, unsigned ncols>
00528 bool
00529 vnl_matrix_fixed<T,nrows,ncols>::is_identity() const
00530 {
00531 T const zero(0);
00532 T const one(1);
00533 for (unsigned int i = 0; i < nrows; ++i)
00534 for (unsigned int j = 0; j < ncols; ++j)
00535 {
00536 T xm = this->data_[i][j];
00537 if ( !((i == j) ? (xm == one) : (xm == zero)) )
00538 return false;
00539 }
00540 return true;
00541 }
00542
00543
00544 template <class T, unsigned nrows, unsigned ncols>
00545 bool
00546 vnl_matrix_fixed<T,nrows,ncols>::is_identity(double tol) const
00547 {
00548 T one(1);
00549 for (unsigned int i = 0; i < nrows; ++i)
00550 for (unsigned int j = 0; j < ncols; ++j)
00551 {
00552 T xm = this->data_[i][j];
00553 abs_t absdev = (i == j) ? vnl_math_abs(xm - one) : vnl_math_abs(xm);
00554 if (absdev > tol)
00555 return false;
00556 }
00557 return true;
00558 }
00559
00560 template <class T, unsigned nrows, unsigned ncols>
00561 bool
00562 vnl_matrix_fixed<T,nrows,ncols>::is_zero() const
00563 {
00564 T const zero(0);
00565 for (unsigned int i = 0; i < nrows; ++i)
00566 for (unsigned int j = 0; j < ncols; ++j)
00567 if ( !( this->data_[i][ j] == zero) )
00568 return false;
00569
00570 return true;
00571 }
00572
00573 template <class T, unsigned nrows, unsigned ncols>
00574 bool
00575 vnl_matrix_fixed<T,nrows,ncols>::is_zero(double tol) const
00576 {
00577 for (unsigned int i = 0; i < nrows; ++i)
00578 for (unsigned int j = 0; j < ncols; ++j)
00579 if (vnl_math_abs(this->data_[i][j]) > tol)
00580 return false;
00581
00582 return true;
00583 }
00584
00585 template <class T, unsigned nrows, unsigned ncols>
00586 bool
00587 vnl_matrix_fixed<T,nrows,ncols>::has_nans() const
00588 {
00589 for (unsigned int i = 0; i < nrows; ++i)
00590 for (unsigned int j = 0; j < ncols; ++j)
00591 if (vnl_math_isnan(this->data_[i][j]))
00592 return true;
00593
00594 return false;
00595 }
00596
00597 template <class T, unsigned nrows, unsigned ncols>
00598 bool
00599 vnl_matrix_fixed<T,nrows,ncols>::is_finite() const
00600 {
00601 for (unsigned int i = 0; i < nrows; ++i)
00602 for (unsigned int j = 0; j < ncols; ++j)
00603 if (!vnl_math_isfinite(this->data_[i][j]))
00604 return false;
00605
00606 return true;
00607 }
00608
00609
00610 template <class T, unsigned nrows, unsigned ncols>
00611 void
00612 vnl_matrix_fixed<T,nrows,ncols>::assert_finite_internal() const
00613 {
00614 if (is_finite())
00615 return;
00616
00617 vcl_cerr << "\n\n" __FILE__ ": " << __LINE__ << ": matrix has non-finite elements\n";
00618
00619 if (rows() <= 20 && cols() <= 20)
00620 vcl_cerr << __FILE__ ": here it is:\n" << *this << '\n';
00621 else
00622 {
00623 vcl_cerr << __FILE__ ": it is quite big (" << rows() << 'x' << cols() << ")\n"
00624 << __FILE__ ": in the following picture '-' means finite and '*' means non-finite:\n";
00625
00626 for (unsigned int i=0; i<rows(); ++i)
00627 {
00628 for (unsigned int j=0; j<cols(); ++j)
00629 vcl_cerr << char(vnl_math_isfinite(this->data_[i][ j]) ? '-' : '*');
00630 vcl_cerr << '\n';
00631 }
00632 }
00633 vcl_cerr << __FILE__ ": calling abort()\n";
00634 vcl_abort();
00635 }
00636
00637
00638 template <class T, unsigned nrows, unsigned ncols>
00639 void
00640 vnl_matrix_fixed<T,nrows,ncols>::assert_size_internal(unsigned rs,unsigned cs) const
00641 {
00642 if (nrows!=rs || ncols!=cs)
00643 {
00644 vcl_cerr << __FILE__ ": size is " << nrows << 'x' << ncols
00645 << ". should be " << rs << 'x' << cs << vcl_endl;
00646 vcl_abort();
00647 }
00648 }
00649
00650 template <class T, unsigned nrows, unsigned ncols>
00651 bool
00652 vnl_matrix_fixed<T,nrows,ncols>::read_ascii(vcl_istream& s)
00653 {
00654 if (!s.good())
00655 {
00656 vcl_cerr << __FILE__ ": vnl_matrix_fixed<T,nrows,ncols>::read_ascii: Called with bad stream\n";
00657 return false;
00658 }
00659
00660 for (unsigned int i = 0; i < nrows; ++i)
00661 for (unsigned int j = 0; j < ncols; ++j)
00662 s >> this->data_[i][j];
00663
00664 return s.good() || s.eof();
00665 }
00666
00667
00668 template <class T, unsigned nrows, unsigned ncols>
00669 vnl_matrix_fixed<T,nrows,ncols>&
00670 vnl_matrix_fixed<T,nrows,ncols>::flipud()
00671 {
00672 for (unsigned int r1 = 0; 2*r1+1 < nrows; ++r1)
00673 {
00674 unsigned int r2 = nrows - 1 - r1;
00675 for (unsigned int c = 0; c < ncols; ++c)
00676 {
00677 T tmp = this->data_[r1][c];
00678 this->data_[r1][c] = this->data_[r2][c];
00679 this->data_[r2][c] = tmp;
00680 }
00681 }
00682 return *this;
00683 }
00684
00685
00686 template <class T, unsigned nrows, unsigned ncols>
00687 vnl_matrix_fixed<T,nrows,ncols>&
00688 vnl_matrix_fixed<T,nrows,ncols>::fliplr()
00689 {
00690 for (unsigned int c1 = 0; 2*c1+1 < ncols; ++c1)
00691 {
00692 unsigned int c2 = ncols - 1 - c1;
00693 for (unsigned int r = 0; r < nrows; ++r)
00694 {
00695 T tmp = this->data_[r][c1];
00696 this->data_[r][c1] = this->data_[r][c2];
00697 this->data_[r][c2] = tmp;
00698 }
00699 }
00700 return *this;
00701 }
00702
00703 template <class T, unsigned nrows, unsigned ncols>
00704 typename vnl_matrix_fixed<T,nrows,ncols>::abs_t
00705 vnl_matrix_fixed<T,nrows,ncols>::operator_one_norm() const
00706 {
00707 abs_t m(0);
00708 for (unsigned int j=0; j<ncols; ++j)
00709 {
00710 abs_t t(0);
00711 for (unsigned int i=0; i<nrows; ++i)
00712 t += vnl_math_abs( this->data_[i][j] );
00713 if (t > m)
00714 m = t;
00715 }
00716 return m;
00717 }
00718
00719 template <class T, unsigned nrows, unsigned ncols>
00720 typename vnl_matrix_fixed<T,nrows,ncols>::abs_t
00721 vnl_matrix_fixed<T,nrows,ncols>::operator_inf_norm() const
00722 {
00723 abs_t m(0);
00724 for (unsigned int i=0; i<nrows; ++i)
00725 {
00726 abs_t t(0);
00727 for (unsigned int j=0; j<ncols; ++j)
00728 t += vnl_math_abs( this->data_[i][j] );
00729 if (t > m)
00730 m = t;
00731 }
00732 return m;
00733 }
00734
00735
00736 template <class T, unsigned nrows, unsigned ncols>
00737 vnl_matrix_fixed<T,nrows,ncols>&
00738 vnl_matrix_fixed<T,nrows,ncols>::inplace_transpose()
00739 {
00740 assert(nrows==ncols);
00741 for (unsigned i = 0; i < nrows; ++i)
00742 for (unsigned j = i+1; j < ncols; ++j)
00743 {
00744 T t = this->data_[i][j];
00745 this->data_[i][j] = this->data_[j][i];
00746 this->data_[j][i] = t;
00747 }
00748 return *this;
00749 }
00750
00751
00752
00753
00754
00755
00756
00757
00758 #if defined(VCL_VC_6) && !defined(__GCCXML__)
00759
00760 template<class VecA, class VecB, class RM>
00761 RM
00762 outer_product_fixed_calc_helper<VecA,VecB,RM>::calc( VecA const& a, VecB const& b )
00763 {
00764 RM out;
00765 for (unsigned int i = 0; i < VecA::SIZE; ++i)
00766 for (unsigned int j = 0; j < VecB::SIZE; ++j)
00767 out[i][j] = a[i] * b[j];
00768 return out;
00769 };
00770
00771 #define VNL_OUTER_PRODUCT_FIXED_INSTANTIATE( T, M, N ) \
00772 template struct outer_product_fixed_calc_helper<vnl_vector_fixed<T,M >, \
00773 vnl_vector_fixed<T,N >, \
00774 vnl_matrix_fixed<T,M,N > >
00775
00776 #else // no need for workaround; declare the function sanely.
00777
00778 template <class T, unsigned m, unsigned n>
00779 vnl_matrix_fixed<T,m,n>
00780 outer_product(vnl_vector_fixed<T,m> const& a, vnl_vector_fixed<T,n> const& b)
00781 {
00782 vnl_matrix_fixed<T,m,n> out;
00783 for (unsigned int i = 0; i < m; ++i)
00784 for (unsigned int j = 0; j < n; ++j)
00785 out[i][j] = a[i] * b[j];
00786 return out;
00787 }
00788
00789 #define VNL_OUTER_PRODUCT_FIXED_INSTANTIATE( T, M, N ) \
00790 template vnl_matrix_fixed<T,M,N > outer_product(vnl_vector_fixed<T,M > const&,\
00791 vnl_vector_fixed<T,N > const& )
00792
00793 #endif // VC60 outer_product workaround
00794
00795
00796 #undef VNL_MATRIX_FIXED_INSTANTIATE
00797 #define VNL_MATRIX_FIXED_INSTANTIATE(T, M, N) \
00798 template class vnl_matrix_fixed<T,M,N >; \
00799 VNL_OUTER_PRODUCT_FIXED_INSTANTIATE( T, M, N )
00800
00801 #endif // vnl_matrix_fixed_txx_