STIR  6.2.0
Filter.h
Go to the documentation of this file.
1 //
2 //
3 
11 /*
12  Copyright (C) 2000 PARAPET partners
13  SPDX-License-Identifier: Apache-2.0
14  Copyright (C) 2000- 2002, IRSL
15  See STIR/LICENSE.txt for details
16 */
17 
18 #ifndef __FILTER_H__
19 #define __FILTER_H__
20 
21 #include "stir/TimedObject.h"
22 #include "stir/Array.h"
23 #include "stir_experimental/fft.h"
24 #include <string>
25 
26 START_NAMESPACE_STIR
27 
28 template <int num_dimensions, typename elemT>
29 class Array;
30 
37 template <class T>
38 class Filter1D : public TimedObject
39 {
40  // AZ public, as need to access length for parallel FBP3DRP code (TODO)
41  // protected:
42 public:
45 
46 public:
47  Filter1D(const int length_v)
48  : filter(1, length_v)
49  {}
50  Filter1D(const Array<1, T>& filter)
51  : filter(filter)
52  {}
53  virtual ~Filter1D() {}
54 
56 
60  inline void apply(Array<1, T>& data);
62  inline void apply(Array<2, T>& data);
64  inline void apply(Array<3, T>& data);
65 
66  virtual std::string parameter_info() const = 0;
67 };
68 
69 // TODO can't be const due to start_timers()
70 template <class T>
71 void
73 {
74  start_timers();
75  const int length = filter.get_length();
76  if (length == 0)
77  return;
78  assert(length > 0);
79 
80  Array<1, T> Padded = data;
81  Padded.set_offset(1);
82  Padded.grow(1, length);
83  convlvC(Padded, filter, length);
84 
85  Padded.set_offset(data.get_min_index());
86 
87  for (int j = data.get_min_index(); j <= data.get_max_index(); j++)
88  data[j] = Padded[j];
89 
90  stop_timers();
91 }
92 
93 template <class T>
94 void
96 {
97  for (int i = data.get_min_index(); i <= data.get_max_index(); i++)
98  apply(data[i]);
99 }
100 
101 template <class T>
102 void
104 {
105  for (int i = data.get_min_index(); i <= data.get_max_index(); i++)
106  apply(data[i]);
107 }
108 
113 template <class T>
114 class Filter2D : public TimedObject
115 {
116 public:
117  Filter2D(int height_v, int width_v)
118  : height(height_v),
119  width(width_v),
120  filter(1, width_v * height_v)
121  {}
122 
123  virtual ~Filter2D() {}
124  inline void apply(Array<1, T>& data) const;
125 
126  // TODO ???
127  void padd_scale_filter(int height_proj, int width_proj);
128 
129  virtual std::string parameter_info() const = 0;
130 
131 protected:
132  int height;
133  int width;
134 
137 };
138 
139 template <class T>
140 void
141 Filter2D<T>::apply(Array<1, T>& data) const
142 {
143 
144  assert(data.get_length() == 2 * height * width);
145 
146  assert(data.get_min_index() == 1);
147 
148  int j, k;
149  /*TODO copy lines from Filter_proj_Colsher for this convertion, and call apply with a Array<2,T>
150  Array<1,T> data(1, 2*height*width);
151 
152  for (int h = 1, j=1; h <= data2d.height(); h++)
153  for (int w=1; w <= data2d.get_width(); w++) {
154  data[j++] = data2d[h][w];
155  }
156  */
157 
158  Array<1, int> nn(1, 2);
159  nn[1] = height;
160  nn[2] = width;
161  fourn(data, nn, 2, 1);
162  for (j = 1, k = 1; j < width * height + 1; j++)
163  {
164  // KT 11/04/2000 normalise with total length
165  // divide by width*height to take Num Recipes convention into account that the
166  // 'inverse' Fourier transform needs scaling by the total number of data points.
167  data[k++] *= filter[j] / (width * height);
168  data[k++] *= filter[j] / (width * height);
169  }
170  fourn(data, nn, 2, -1);
171 };
172 
173 END_NAMESPACE_STIR
174 
175 #endif
int get_length() const
return number of elements in this vector
Definition: VectorWithOffset.inl:534
Array< 1, float > filter
Stores the filter in frequency space (will be private sometime)
Definition: Filter.h:44
int get_min_index() const
get value of first valid index
Definition: VectorWithOffset.inl:116
void fourn(Array< 1, float > &data, Array< 1, int > &nn, int ndim, int isign)
n-dimensional FFT
virtual void grow(const IndexRange< num_dimensions > &range)
alias for resize()
Definition: Array.inl:83
2-dimensional filters (filtering done by FFTs)
Definition: Filter.h:114
void convlvC(Array< 1, float > &data, const Array< 1, float > &filter, int n)
Convolve data with a filter which is given in frequency space.
defines the Array class for multi-dimensional (numeric) arrays
base class for all objects which need timers. At the moment, there&#39;s only a CPU timer.
Definition: TimedObject.h:35
Array< 1, T > filter
Stores the filter in the &#39;funny&#39; but efficient Numerical Recipes format.
Definition: Filter.h:136
declares the stir::TimedObject class
Preliminary class for 1D filtering using FFTs.
Definition: Filter.h:38
int get_max_index() const
get value of last valid index
Definition: VectorWithOffset.inl:123
void set_offset(const int min_index)
change value of starting index
Definition: VectorWithOffset.inl:332
Declaration of FFT routines.