sparse_image_solver.hh
Go to the documentation of this file.
1 /* -*- mia-c++ -*-
2  *
3  * This file is part of MIA - a toolbox for medical image analysis
4  * Copyright (c) Leipzig, Madrid 1999-2013 Gert Wollny
5  *
6  * MIA is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 3 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with MIA; if not, see <http://www.gnu.org/licenses/>.
18  *
19  */
20 
21 #ifndef mia_2d_sparse_image_solver_hh
22 #define mia_2d_sparse_image_solver_hh
23 
24 #include <mia/core/factory.hh>
26 #include <mia/2d/image.hh>
27 
29 
32 
34 typedef std::shared_ptr<C2DImageSparseSolver > P2DImageSparseSolver;
35 
38 
45 struct EXPORT_2D C2DImageSparseSolverTestPath {
46  C2DImageSparseSolverTestPath();
47 };
48 
50 
53 
55 typedef std::shared_ptr<C2DImageSolverAmultx> P2DImageSolverAmultx;
56 
59 
62 
63 
64 
74 template <typename T>
75 struct multiply<T2DImage<T> > {
82  static void apply(T2DImage<T>& result, const typename TSparseSolver<T2DImage<T> >::A_mult_x& A, const T2DImage<T>& x);
83 };
84 
85 template <typename T>
86 void multiply<T2DImage<T> >::apply(T2DImage<T>& result, const typename TSparseSolver<T2DImage<T> >::A_mult_x& A, const T2DImage<T>& X)
87 {
88  assert(result.get_size() == X.get_size());
89  assert(result.get_size() == A.get_size());
90 
91  long b = A.get_boundary_size();
92  long nx = X.get_size().x - 2 * b;
93  long ny = X.get_size().y - 2 * b;
94  copy(X.begin(), X.begin() + b * X.get_size().x + b, result.begin());
95  auto ir = result.begin() + b * X.get_size().x + b;
96  auto ix = X.begin() + b * X.get_size().x + b;
97 
98  for (int y = 0; y < ny; ++y) {
99  int x = 0;
100  for (; x < nx; ++x, ++ix, ++ir)
101  *ir = static_cast<T>(A(ix));
102  for (; x < (int)X.get_size().x; ++x, ++ix, ++ir)
103  *ir = *ix;
104  }
105  copy(ix, X.end(), ir);
106 }
107 
112 
113 
115 
116 #endif