/**
  \file

  \brief Graybody implementation.

  \who \kpr
*/
#include "Detector.h"
#include "Graybody.h"

namespace mu = mutils;

namespace manticore {

/// Detector support.
namespace graybody {

extern "C" {

/// Parameters required by flux().
struct FluxData {
  double T,      ///< Temperature.
         Sigma;  ///< Mass surface density.

  const Detector *d;  ///< Detector instance (may be `nullptr`).

  const Graybody *g;  ///< Graybody instance.

  double (Graybody::*I)(double, double, double) const;
};

/// Flux density integrand.
/// \param[in] nu  Frequency (Hz).
/// \param[in] fd  FluxData instance.
double flux(double nu, void *fd) {
  auto &p = *static_cast<const FluxData*>(fd);
  return (p.g->*(p.I))(nu, p.T, p.Sigma)*(p.d ? p.d->response(nu) : 1.0);
}

}  // extern "C"

}  // namescape graybody


double Graybody::Itot(double (Graybody::*I)(double, double, double) const,
                      double T, double Sigma, const Detector *detect,
                      double sr) const
{
  std::pair<double, double> freqRange =
    (detect ? detect->freqRange() :
              std::pair<double,double>{0.0,12.0*k_Boltzman*T/h_Planck});

  double result, abserr;
  graybody::FluxData fd = {T, Sigma, detect, this, I};
  gsl_function f = {graybody::flux, &fd};
  gsl_integration_qag(&f, freqRange.first, freqRange.second,
                      0.0, err_, subLimit, intKey,
                      integ_, &result, &abserr);
  return result*sr;
}

}  // namespace manticore