Cheetah - SKA - PSS - Prototype Time Domain Search Pipeline
ConvolvePlan.cpp
1 /*
2  * The MIT License (MIT)
3  *
4  * Copyright (c) 2016 The SKA organisation
5  *
6  * Permission is hereby granted, free of charge, to any person obtaining a copy
7  * of this software and associated documentation files (the "Software"), to deal
8  * in the Software without restriction, including without limitation the rights
9  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
10  * copies of the Software, and to permit persons to whom the Software is
11  * furnished to do so, subject to the following conditions:
12  *
13  * The above copyright notice and this permission notice shall be included in all
14  * copies or substantial portions of the Software.
15  *
16  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
17  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
18  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
19  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
20  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
21  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
22  * SOFTWARE.
23  */
24 #include "cheetah/utils/ConvolvePlan.h"
25 
26 
27 namespace ska {
28 namespace cheetah {
29 namespace utils {
30 
31 
32 ConvolvePlan::ConvolvePlan(std::size_t npts, float* a, float* b, float* output)
33  : _a(a)
34  , _b(b)
35  , _output(output)
36  , _npts(npts)
37  , _c(fftwf_alloc_complex(npts))
38 {
39  if(_c == nullptr) throw std::bad_alloc();
40  _d=fftwf_alloc_complex(npts);
41  if(_d == nullptr) {
42  fftwf_free(_c);
43  throw std::bad_alloc();
44  }
45  _fwda = fftwf_plan_dft_r2c_1d(npts, _a, _c, FFTW_MEASURE);
46  if(_fwda == nullptr) {
47  fftwf_free(_c);
48  fftwf_free(_d);
49  throw std::bad_alloc();
50  }
51  _fwdb = fftwf_plan_dft_r2c_1d(npts, _b, _d, FFTW_MEASURE);
52  if(_fwdb == nullptr) {
53  fftwf_free(_c);
54  fftwf_free(_d);
55  fftwf_destroy_plan(_fwda);
56  throw std::bad_alloc();
57  }
58  _bwd = fftwf_plan_dft_c2r_1d(npts, _c, _output, FFTW_MEASURE);
59  if(_bwd == nullptr) {
60  fftwf_free(_c);
61  fftwf_free(_d);
62  fftwf_destroy_plan(_fwda);
63  fftwf_destroy_plan(_fwdb);
64  throw std::bad_alloc();
65  }
66  _scale = (float) npts;
67 }
68 
69 ConvolvePlan::~ConvolvePlan()
70 {
71  fftwf_free(_c);
72  fftwf_free(_d);
73  fftwf_destroy_plan(_fwda);
74  fftwf_destroy_plan(_fwdb);
75  fftwf_destroy_plan(_bwd);
76 }
77 
79 {
80  std::size_t i, j;
81  fftwf_execute(_fwda);
82  fftwf_execute(_fwdb);
83  for(i = 0; i < _npts; ++i) {
84  for(j = 0; j < 2; ++j) {
85  _c[i][j] *= _d[i][j]/_scale;
86  }
87  }
88  fftwf_execute(_bwd);
89 }
90 
91 } // namespace utils
92 } // namespace cheetah
93 } // namespace ska
Some limits and constants for FLDO.
Definition: Brdz.h:35
ConvolvePlan(std::size_t npts, float *a, float *b, float *output)
fast convolution of two data sets a and b of length npts output = a * b