24 #include "../RfiGenerator.h" 26 #pragma GCC diagnostic push 27 #pragma GCC diagnostic ignored "-Wunused-parameter" 28 #include <boost/accumulators/accumulators.hpp> 29 #include <boost/accumulators/statistics/stats.hpp> 30 #include <boost/accumulators/statistics/mean.hpp> 31 #include <boost/accumulators/statistics/variance.hpp> 32 #pragma GCC diagnostic pop 36 namespace generators {
40 RfiGenerator<T>::RfiGenerator()
41 : TimeFrequencyGenerator<T>()
46 RfiGenerator<T>::~RfiGenerator()
51 void RfiGenerator<T>::next(DataType& data)
53 FlaggedDataType flagged_data(data);
54 this->
next(flagged_data);
58 void RfiGenerator<T>::rfi_gaussian_block(FlaggedDataType& flagged_data
59 , std::size_t min_channel
60 , std::size_t max_channel
61 , std::size_t min_sample_number
62 , std::size_t max_sample_number
66 DataType& data =
static_cast<DataType&
>(flagged_data);
67 boost::accumulators::accumulator_set<double, boost::accumulators::stats<boost::accumulators::tag::variance,boost::accumulators::tag::mean> > acc;
68 for_each(data.begin(), data.end(), std::bind<void>(std::ref(acc), std::placeholders::_1));
70 float mean = boost::accumulators::mean(acc);
71 float variance = boost::accumulators::variance(acc);
72 float stddev = sqrtf(variance);
78 float max = mean + sigma * stddev;
81 int centre_sample = int((max_sample_number + min_sample_number)/2);
82 int centre_channel = int((max_channel + min_channel)/2);
83 float channel_gaussian_width = (max_channel - min_channel)/2.0;
84 float time_gaussian_width = ((max_sample_number - min_sample_number)/2.0);
85 float channel_gaussian_variance = powf(channel_gaussian_width/2.355, 2.0f);
86 float time_gaussian_variance = powf(time_gaussian_width/2.355, 2.0f);
87 auto & flags = flagged_data.rfi_flags();
89 for (
unsigned time_sample = min_sample_number; time_sample < max_sample_number; ++time_sample)
91 auto spectrum = data.spectrum(time_sample);
92 auto channel_it = spectrum.begin();
93 channel_it += min_channel;
94 unsigned channel_number = min_channel;
95 auto flag_slice = flags.spectrum(time_sample);
96 auto flag_it = flag_slice.begin();
97 flag_it += min_channel;
98 while(channel_number <= max_channel)
100 *channel_it += max * std::exp(-( (powf(
float(channel_number) -
float(centre_channel), 2.0f)/(2.0f*channel_gaussian_variance))
101 + (powf(
float(time_sample) -
float(centre_sample), 2.0f)/(2.0f*time_gaussian_variance))));
102 if(time_sample<centre_sample+std::floor(time_gaussian_width)+1 &&
103 time_sample>centre_sample-std::floor(time_gaussian_width)-1 &&
104 channel_number<centre_channel+std::floor(channel_gaussian_width)+1 &&
105 channel_number>centre_channel-std::floor(channel_gaussian_width)-1)
118 void RfiGenerator<T>::rfi_ramp_block(FlaggedDataType& flagged_data
119 , std::size_t min_channel
120 , std::size_t max_channel
121 , std::size_t min_sample_number
122 , std::size_t max_sample_number)
124 DataType& data =
static_cast<DataType&
>(flagged_data);
125 boost::accumulators::accumulator_set<double, boost::accumulators::stats<boost::accumulators::tag::variance,boost::accumulators::tag::mean> > acc;
126 for_each(data.begin(), data.end(), std::bind<void>(std::ref(acc), std::placeholders::_1));
128 float mean = boost::accumulators::mean(acc);
129 float variance = boost::accumulators::variance(acc);
130 float stddev = sqrtf(variance);
135 float max = mean + 10.0f * stddev;
136 float slope = 1.0f/float(max_channel - min_channel);
138 auto & flags = flagged_data.rfi_flags();
140 for (
unsigned time_sample = min_sample_number; time_sample < max_sample_number; ++time_sample)
142 auto spectrum = data.spectrum(time_sample);
143 auto channel_it = spectrum.begin();
144 channel_it += min_channel;
145 unsigned channel_number = min_channel;
146 auto flag_slice = flags.spectrum(time_sample);
147 auto flag_it = flag_slice.begin();
148 flag_it += min_channel;
149 while(channel_number <= max_channel)
151 *channel_it += max * slope * (channel_number - min_channel);
161 void RfiGenerator<T>::rfi_broadband_block(FlaggedDataType& flagged_data
162 , std::size_t min_sample_number
163 , std::size_t max_sample_number
167 DataType& data =
static_cast<DataType&
>(flagged_data);
168 boost::accumulators::accumulator_set<double, boost::accumulators::stats<boost::accumulators::tag::variance,boost::accumulators::tag::mean> > acc;
169 for_each(data.begin(), data.end(), std::bind<void>(std::ref(acc), std::placeholders::_1));
171 float mean = boost::accumulators::mean(acc);
172 float variance = boost::accumulators::variance(acc);
173 float stddev = sqrtf(variance);
174 std::size_t min_channel = 1;
175 std::size_t max_channel = flagged_data.number_of_channels()-1;
178 float width = float((max_sample_number-min_sample_number)/2.0f);
179 float variance_of_gaussian = powf(width/2.355, 2.0f);
185 float max = mean + sigma * stddev;
188 int centre_sample = int((max_sample_number + min_sample_number)/2);
190 auto & flags = flagged_data.rfi_flags();
192 for (
unsigned time_sample = min_sample_number; time_sample < max_sample_number; ++time_sample)
194 auto spectrum = data.spectrum(time_sample);
195 auto channel_it = spectrum.begin();
196 channel_it += min_channel;
197 unsigned channel_number = min_channel;
198 auto flag_slice = flags.spectrum(time_sample);
199 auto flag_it = flag_slice.begin();
200 flag_it += min_channel;
201 while(channel_number <= max_channel)
203 *channel_it += max * std::exp(-(powf(
float(time_sample) -
float(centre_sample), 2.0f))/(2.0f*variance_of_gaussian));
204 if(time_sample<centre_sample+std::floor(width)+1 &&
205 time_sample>centre_sample-std::floor(width)-1)
virtual void next(DataType &data) override
do nothing with the passed data
Some limits and constants for FLDO.