24 #include "cheetah/data/SpCcl.h" 25 #include "cheetah/utils/ModifiedJulianClock.h" 26 #include "pss/astrotypes/units/Time.h" 34 static const boost::units::quantity<MilliSeconds, double> dm_const_factor = (4.15 * 1000000.0 * boost::units::si::milli * boost::units::si::seconds);
36 template<
typename NumericalRep>
37 SpCcl<NumericalRep>::SpCcl()
38 : _dm_range(Dm(
std::numeric_limits<typename Dm::value_type>::max() * parsecs_per_cube_cm),
39 Dm(
std::numeric_limits<typename Dm::value_type>::min() * parsecs_per_cube_cm))
43 template<
typename NumericalRep>
44 SpCcl<NumericalRep>::SpCcl(SpCcl<NumericalRep>::BlocksType
const& blocks, data::DimensionIndex<data::Time> offset)
45 : _dm_range(Dm(
std::numeric_limits<typename Dm::value_type>::max() * parsecs_per_cube_cm),
46 Dm(
std::numeric_limits<typename Dm::value_type>::min() * parsecs_per_cube_cm))
50 auto const& first_block = _blocks.front();
51 _offset_time = pss::astrotypes::units::duration_cast<boost::units::quantity<data::MilliSeconds,double>>((double)offset * first_block->sample_interval());
52 _start_time=first_block->start_time() + _offset_time;
53 std::pair<typename TimeFrequencyType::FrequencyType, typename TimeFrequencyType::FrequencyType> freq = first_block->low_high_frequencies();
54 _dm_delay_factor = ((1.0/(freq.first * freq.first)) - (1.0/(freq.second * freq.second)));
58 template<
typename NumericalRep>
59 SpCcl<NumericalRep>::~SpCcl()
63 template<
typename NumericalRep>
69 template<
typename NumericalRep>
75 template<
typename NumericalRep>
78 if(dm > _dm_range.second) _dm_range.second =
dm;
79 if(dm < _dm_range.first) _dm_range.first =
dm;
82 template<
typename NumericalRep>
85 dm_max_min(cand.
dm());
87 auto it = std::find_if(BaseT::rbegin()
91 return start_time >= cand.
tstart();
93 auto r_it = BaseT::insert(it.base(), cand);
94 r_it->tend(start_time + dm_const_factor * (cand.
dm() * _dm_delay_factor).value());
97 template<
typename NumericalRep>
100 dm_max_min(cand.
dm());
102 auto it = std::find_if(BaseT::rbegin()
106 return start_time >= cand.
tstart();
108 BaseT::insert(it.base(), cand);
111 template<
typename NumericalRep>
114 dm_max_min(cand.dm());
115 if(BaseT::empty() || cand.tstart() >= BaseT::back().tstart()) {
116 BaseT::emplace_back(std::move(cand));
120 auto it = std::lower_bound(BaseT::begin()
125 return start_time > cand.
tstart();
127 BaseT::insert(it,std::move(cand));
131 template<
typename NumericalRep>
134 dm_max_min(cand.dm());
136 auto it = std::lower_bound(BaseT::begin()
139 , [&](
SpCandidateType const& cand, boost::units::quantity<MilliSeconds, double> start_time)
141 return start_time > cand.
tstart();
144 BaseT::insert(it, std::move(cand));
147 template<
typename NumericalRep>
150 dm_max_min(cand.dm());
152 cand.tend(start_time + dm_const_factor * (_dm_delay_factor * cand.dm()).value());
153 auto it = std::lower_bound(BaseT::begin()
156 , [&](
SpCandidateType const& cand, boost::units::quantity<MilliSeconds, double> start_time)
158 return start_time > cand.
tstart();
160 BaseT::insert(it,std::move(cand));
163 template<
typename NumericalRep>
165 boost::units::quantity<MilliSeconds, double>
start_time,
166 boost::units::quantity<MilliSeconds, double> pulse_width,
171 auto it = std::lower_bound(BaseT::begin()
174 , [&](
SpCandidateType const& cand, boost::units::quantity<MilliSeconds, double> start_time)
176 return start_time > cand.
tstart();
179 ,start_time + dm_const_factor * (_dm_delay_factor * dm).value()
183 template<
typename NumericalRep>
185 boost::units::quantity<MilliSeconds, double>
start_time,
186 boost::units::quantity<MilliSeconds, double> pulse_width,
191 BaseT::emplace_back(dm, start_time, pulse_width
192 ,start_time + dm_const_factor * (_dm_delay_factor * dm).value()
196 template<
typename NumericalRep>
199 if (BaseT::begin() == BaseT::end())
205 return ConstDataIterator(window, BaseT::begin(), BaseT::end(), _blocks, _start_time);
209 template<
typename NumericalRep>
215 template<
typename NumericalRep>
221 template<
typename NumericalRep>
224 return _start_time + pss::astrotypes::units::duration_cast<utils::julian_day>(cand.
tstart());
227 template<
typename NumericalRep>
233 template<
typename NumericalRep>
234 template<
typename PredicateT>
237 Dm low = _dm_range.second;
238 Dm high = _dm_range.first;
239 auto it = std::remove_if(begin(), end(),
243 if(cand.
dm() > high) high = cand.
dm();
244 if(cand.
dm() < low) low = cand.
dm();
249 if (it != BaseT::end())
251 erase(it, BaseT::end());
253 _dm_range.first = low;
254 _dm_range.second = high;
258 template<
typename NumericalRep>
260 , ConstCandidateIterator cand_it_begin
261 , ConstCandidateIterator
const& cand_it_end
262 , SpCcl::BlocksType
const& tf_blocks
263 , utils::ModifiedJulianClock::time_point
const& epoch
266 , _cand_it(cand_it_begin)
267 , _cand_it_end(cand_it_end)
268 , _tf_it(tf_blocks.begin())
269 , _tf_end(tf_blocks.end())
270 , _tf_start_it(tf_blocks.begin())
273 if (_cand_it != _cand_it_end)
275 _span = ((*_tf_it)->number_of_spectra() *
static_cast<boost::units::quantity<MilliSeconds,double>
>((*_tf_it)->sample_interval()));
277 this->next_candidate();
281 template<
typename NumericalRep>
284 if( _cand_it == _cand_it_end )
288 auto increment_block = [
this]()
291 if (_tf_it == _tf_end)
295 _span = ((*_tf_it)->number_of_spectra() *
static_cast<boost::units::quantity<MilliSeconds,double>
>((*_tf_it)->sample_interval()));
300 boost::units::quantity<MilliSeconds, double>
start_time = _cand_it->tstart() - _window.ms_before();
302 if (start_time < _duration - _span)
304 start_time = _duration - _span;
307 auto end_time = _cand_it->tend() + _window.ms_after();
308 auto start_cand = _cand_it;
311 while(++_cand_it != _cand_it_end)
313 if(_cand_it->tstart() - _window.ms_before() <= end_time) {
314 auto next_end_time = _cand_it->tend() + _window.ms_after();
315 if(next_end_time > end_time) end_time = next_end_time;
322 _cands.emplace_back(start_cand, _cand_it);
324 while(start_time >= _duration)
326 if (!increment_block())
328 throw panda::Error(
"Candidate Start time larger than end of data!!");
332 _start_time = _epoch + pss::astrotypes::units::duration_cast<utils::julian_day>(
start_time);
334 data::DimensionIndex<data::Time> start_offset((*_tf_it)->number_of_spectra() - (std::size_t)(((_duration - start_time)/
static_cast<boost::units::quantity<MilliSeconds,double>
>((*_tf_it)->sample_interval()))));
336 PANDA_LOG_DEBUG <<
"Start offset:" << start_offset;
337 while (end_time >= _duration)
339 data::DimensionIndex<data::Time> end_offset((std::size_t)((*_tf_it)->number_of_spectra()-1));
340 PANDA_LOG_DEBUG <<
"End Offset:" << end_offset;
341 _slices.push_back((*_tf_it)->slice(data::DimensionSpan<data::Time>(start_offset, end_offset)));
342 start_offset = data::DimensionIndex<data::Time>(0);
343 if (!increment_block())
345 _slices.back()._last_block =
true;
349 data::DimensionIndex<data::Time> end_offset((*_tf_it)->number_of_spectra() - 1 - (std::size_t)(((_duration - end_time)/
static_cast<boost::units::quantity<MilliSeconds,double>
>((*_tf_it)->sample_interval()))));
350 PANDA_LOG_DEBUG <<
"Start Offset=" << start_offset <<
" End Offset:" << end_offset;
351 _slices.emplace_back((*_tf_it)->slice(data::DimensionSpan<data::Time>(start_offset,end_offset)),
true);
355 template<
typename NumericalRep>
358 return _slices.front().last_block();
361 template<
typename NumericalRep>
367 template<
typename NumericalRep>
370 auto& slice = _slices.front();
371 _start_time += pss::astrotypes::units::duration_cast<utils::julian_day>(slice.slice().template dimension<data::Time>() * (*_tf_start_it)->sample_interval());
373 if(slice.last_block()) {
378 if(_slices.size() == 0)
385 template<
typename NumericalRep>
388 return _cands.front().first;
391 template<
typename NumericalRep>
394 return _cands.
front().second;
397 template<
typename NumericalRep>
405 template<
typename NumericalRep>
408 return _slices.
front().slice();
411 template<
typename NumericalRep>
414 if(d._slices.size() != this->_slices.size())
419 if(d._slices.size() == 0)
return true;
421 return d._slices.front().slice() == this->_slices.front().slice();
424 template<
typename NumericalRep>
427 return !(d == *
this);
ValueType const & front() const
the first emelment
ConstDataIterator data_begin(CandidateWindow const &window) const
Extract data corresponding to the candidates start and end times.
A simple record to hold 'candidate' proprerties.
void push_back(SpCandidateType const &cand)
Push back method for inserting Single pulse Candidates but without calculating tend They will be ins...
void push_back_calculate_duration(SpCandidateType const &cand)
Push back method for inserting Single pulse Candidates calculating temd based on dispersion and freq ...
std::pair< Dm, Dm > const & dm_range() const
the max and minimum Dispersion Measure values of the candidates
BlocksType const & tf_blocks() const
Access the TimeFrequency blocks associated with the single pulse candidtes.
void emplace(SpCandidateType &&cand)
Emplace back method for inserting Single pulse Candidates, ensuring ordering by start_time. They will be inserted in a sorted order with respect to the start time of each candidate instead of insert to the end.
TimePointType const & start_time() const
Get the start time of the candidate.
bool end_block() const
retunrs true is the current block is the end of contiguous series
void emplace_back(SpCandidateType &&cand)
Add a Single pulse Candidate to the end of the cand list, without calculating the duration Checking ...
Dm const & dm() const
access a reference to dm.
Some limits and constants for FLDO.
ConstCandidateIterator candidate_begin() const
retunrs true is the current block is the end of contiguous series
MsecTimeType const & tstart() const
Get the start time of the candidate.
utils::ModifiedJulianClock::time_point const & start_time() const
the start time of the first relevant data (taking into account the offset)
ConstDataIterator & operator++()
pre-increment operator
Dm const & dm() const
access a reference to dm.
define a window of data surrounding a single pulse candidate
boost::units::quantity< data::MilliSeconds, double > offset_time() const
the time duration from the beginning of the first block to start_time()
Sigma const & sigma() const
access a reference to standard deviation (sigma).
void remove_if(PredicateT const &)
remove candidates that meet the specified condition
utils::ModifiedJulianClock::time_point const & start_time() const
: Get the start time