PowerSpectrum.h

#ifndef _POWER_SPECTRUM_H_
#define _POWER_SPECTRUM_H_
#include <vector>
#include "mspass/utility/Metadata.h"
#include "mspass/utility/MsPASSError.h"
#include "mspass/utility/ErrorLogger.h"
#include "mspass/seismic/BasicSpectrum.h"
namespace mspass::seismic
{
class PowerSpectrum : public mspass::seismic::BasicSpectrum,
                      public mspass::utility::Metadata
{
public:
  std::string spectrum_type;
  std::vector<double> spectrum;
  mspass::utility::ErrorLogger elog;
  PowerSpectrum();
  template <class T> PowerSpectrum(const mspass::utility::Metadata& md,
    const std::vector<T>& d,
      const double dfin,
        const std::string nm,
          const double f0in,
            const double dtin,
              const int npts_in);
  PowerSpectrum(const PowerSpectrum& parent);
  PowerSpectrum& operator=(const PowerSpectrum& parent);
  PowerSpectrum& operator+=(const PowerSpectrum& other);
  std::vector<double> amplitude() const;
  double power(const double f) const;
 
  double frequency(const int sample_number) const
  {
    const std::string base_error("PowerSpectrum::frequency:  ");
    if(sample_number<0) throw mspass::utility::MsPASSError(base_error
        + "Sample number parameter passed cannot be negative");
    if(sample_number>=(this->nf())) throw mspass::utility::MsPASSError(base_error
        + "Sample number parameter passed xceeds range of spectrum array");
    return this->f0()+sample_number*this->df();
  };
  std::vector<double> frequencies() const;
  /* \brief Return number of frequencies in estimate (size of spectrum vector).
 
  Users should be aware that when zero padding is used the size of the
  spectrum vector will be larger than half the number of time series samples.*/
  size_t nf()const{return spectrum.size();};
  double Nyquist() const {
    return 1.0/(2.0*this->parent_dt);
  };
};
template <class T> PowerSpectrum::PowerSpectrum(const mspass::utility::Metadata& md,
    const std::vector<T>& d,const double dfin,const std::string nm,
      const double f0in, const double dtin, const int npts_in)
      : BasicSpectrum(dfin,f0in,dtin,npts_in),
          mspass::utility::Metadata(md),
            elog()
{
  spectrum_type=nm;
  spectrum.reserve(d.size());
  for(size_t k=0;k<d.size();++k)
    spectrum.push_back(static_cast<double>(d[k]));
  this->set_live();
};
}  //end namespace
#endif