9#include "NumericalAtmosApprox.h"
10#include "NumericalApproximation.h"
11#include "IException.h"
59 double trap[maxits+1];
70 for(
int i = 0; i < maxits; i++) {
75 interp.
AddData(5, &h[i-4], &trap[i-4]);
80 if(fabs(dss) <= epsilon * fabs(ss))
return ss;
81 if(fabs(dss) <= epsilon2)
return ss;
95 "NumericalAtmosApprox::RombergsMethod() - Caught the following error: ",
99 "NumericalAtmosApprox::RombergsMethod() - Failed to converge in "
100 +
IString(maxits) +
" iterations.",
152 return (0.5 * (b - a) * (begin + end));
156 double delta, tnm, x, sum;
157 it = (int)(pow(2.0, (
double)(n - 2)));
159 delta = (b - a) / tnm;
162 for(
int i = 0; i < it; i++) {
171 return (0.5 * (s + (b - a) * sum / tnm));
177 "NumericalAtmosApprox::RefineExtendedTrap() - Caught the following error: ",
206 am->p_atmosPhi = phi;
207 am->p_atmosCosphi = cos((
PI / 180.0) * phi);
217 "NumericalAtmosApprox::OutrFunc2Bint() - Caught the following error: ",
264 ema = acos(mu) * (180.0 /
PI);
265 sine = sin(ema * (
PI / 180.0));
266 if((am->p_atmosAtmSwitch == 0) || (am->p_atmosAtmSwitch == 2)) {
267 alpha = am->p_atmosSini * sine * am->p_atmosCosphi +
268 am->p_atmosMunot * mu;
271 alpha = am->p_atmosSini * sine * am->p_atmosCosphi -
272 am->p_atmosMunot * mu;
274 phase = acos(alpha) * (180.0 /
PI);
277 if(am->p_atmosAtmSwitch == 0) {
279 result = mu * am->p_atmosPM->
CalcSurfAlbedo(phase, am->p_atmosInc, ema);
283 xx = -am->
AtmosTau() / max(am->p_atmosMunot, 1.0e-30);
293 xx = -am->
AtmosTau() / max(mu, 1.0e-30);
304 if(mu == am->p_atmosMunot) {
305 tfac = am->
AtmosTau() * emunot / (am->p_atmosMunot * am->p_atmosMunot);
308 tfac = (emunot - emu) / (am->p_atmosMunot - mu);
310 if(am->p_atmosAtmSwitch == 1) {
311 result = mu * (phasefn - 1.0) * tfac;
313 else if(am->p_atmosAtmSwitch == 2) {
314 result = am->p_atmosMunot * mu * (phasefn - 1.0) * (1.0 - emunot * emu) / (am->p_atmosMunot + mu);
316 else if(am->p_atmosAtmSwitch == 3) {
317 result = -sine * am->p_atmosCosphi * (phasefn - 1.0) * tfac;
320 string msg =
"NumericalAtmosApprox::InrFunc2Bint() - Invalid value of atmospheric ";
321 msg +=
"switch used as argument to this function";
Isotropic atmos scattering model.
double AtmosTau() const
Return atmospheric Tau value.
double AtmosHga() const
Return atmospheric Hga value.
@ Programmer
This error is for when a programmer made an API call that was illegal.
ErrorType errorType() const
Returns the source of the error for this exception.
Adds specific functionality to C++ strings.
double Evaluate(const double a, const ExtrapType &etype=ThrowError)
Calculates interpolated or extrapolated value of tabulated data set for given domain value.
@ PolynomialNeville
Polynomial interpolation using Neville's algorithm.
void Reset()
Resets the state of the object.
NumericalApproximation(const NumericalApproximation::InterpType &itype=CubicNatural)
Default constructor creates NumericalApproximation object.
@ Extrapolate
Evaluate() attempts to extrapolate if a is outside of the domain. This is only valid for NumericalApp...
vector< double > PolynomialNevilleErrorEstimate()
Retrieves the error estimate for the Neville's polynomial interpolation type.
void AddData(const double x, const double y)
Add a datapoint to the set.
static double InrFunc2Bint(AtmosModel *am, double mu)
Inner function to be integrated.
double RombergsMethod(AtmosModel *am, IntegFunc sub, double a, double b)
This variation on the NumericalApproximation method integrates a specified AtmosModel function rather...
NumericalAtmosApprox(const NumericalApproximation::InterpType &itype=CubicNatural)
Uses Isis::NumericalApproximation constructor.
static double OutrFunc2Bint(AtmosModel *am, double phi)
This function is the outer integrand over mu at specified phi.
double RefineExtendedTrap(AtmosModel *am, IntegFunc sub, double a, double b, double s, unsigned int n)
This variation on the NumericalApproximation method integrates a specified AtmosModel function rather...
IntegFunc
This enum defines function to be integrated by Romberg's method.
@ InnerFunction
Indicates that Romberg's method will integrate the function InrFunc2Bint()
double CalcSurfAlbedo(double pha, double inc, double ema)
Calculate the surface brightness using photometric angle information.
This is free and unencumbered software released into the public domain.
const double PI
The mathematical constant PI.
Namespace for the standard library.