astro
Class OrbitModel

java.lang.Object
  vec_math.AbstractDataModel
      vec_math.AbstractGradientModel
          astro.OrbitModel

All Implemented Interfaces:

DataModel, GradientModel

public class OrbitModel
extends AbstractGradientModel

The orbit model is the differentiable model that fits to a radial velocity of an SB1 orbit. It can be used in a Levenberg-Marquardt schema to solve a Keplerian orbit. The model that fits an SB-1 orbit is given by

       v_r=v_0+e·K*cos(ω)+K*cos(ω+ν),
       

where v_0 is the baricentric velocity of the system, K is the radial velocity amplitude, e is the excentricity of the system, ω is the argument of perihelion, and ν is the true anomaly. From the eccentricity and the eccentric anomaly, the true anomaly follows as

       tan(ν/2) = sqrt((1+e)/(1-e))·tan(E/2)
       

From the time of perihelion passage, the orbital period, and the numeric eccentricity from the orbit, get the excentric anomaly for the current time. Solve Kepler's equation

       M=2π/P·(t-t0)=E-e·sin(E)
       

which is done iteratively.

For calculating the quadrature of an orbit, saee the special class.

Nested Class Summary

static class	OrbitModel.Calculate Supply the orbit parameters t0, P, e, omega, K1 and vB to the command line to get the orbit solution for the specified time.
static class	OrbitModel.Quadrature Calculates the phases of quadrature, defined as ν=±90, or, if a second argument is given, which is then the argument of periastron, as the times of maximum rv, which follows from the defintion of rv to
static class	OrbitModel.Test For testing.

Nested classes/interfaces inherited from class vec_math.AbstractDataModel

AbstractDataModel.Construct

Field Summary

private static double	ACCURACY Accuracy to reach in Kepler's equation.
private double	erradd
private static int	ITERMAX Maximum number of iterations, log_2(2/ACCURACY).

Fields inherited from class vec_math.AbstractDataModel

sigma, x, y

Constructor Summary

OrbitModel(VectorG[] times, double[] y, double[] err, double off)
Measures at differnt point-in-time.

Method Summary

double	evaluateModel(VectorG a, VectorG t) The model that fits an SB-1 orbit is given by
private static double	getExcentricAnomaly(double t, double t0, double p, double e) From the time of perihelion passage, the orbital period, and the numeric eccentricity from the orbit, get the excentric anomaly for the current time.
Matrix	getGradientMatrix(VectorG a) This is the derivative matrix of the orbital solution with respect to the parameters.
double[]	getMeasurementErrors() If we had an external error, we subtract it again.
double[]	getModel(VectorG a) The model that fits an SB-1 orbit is given by
private static double	getOrbitSolution(double t, VectorG a) Calculates the orbit rv for the given time and the specified parameter set.
int	getParameterCount() An orbit has six parameters, namely at indices: 0: The perihelium time, in the units of our time, like JD. 1: The orbital period, in our time units (e.g.
private static double	getTrueAnomaly(double e, double exan) From the eccentricity and the eccentric anomaly return the true anomaly.

Methods inherited from class vec_math.AbstractGradientModel

check, check, getNegativeChi2HalfGradient, getNegativeChi2HalfGradient, getPseudoHessian, getPseudoHessian, LevenbergMarquardtSolver

Methods inherited from class vec_math.AbstractDataModel

getChiSquareModel, getChiSquareModel, getLorentzianModel, getMeasureCount, getMeasures, getModel, getResiduals, getResiduals, getRms, getRobustModel, getTimes, getTotalErrors, LorentzianSolver, RobustSolver, SimplexSolver

Methods inherited from class java.lang.Object

clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait

Methods inherited from interface vec_math.DataModel

getMeasureCount, getMeasures, getResiduals, getTimes, getTotalErrors

Field Detail

ACCURACY

private static double ACCURACY

Accuracy to reach in Kepler's equation.

ITERMAX

private static int ITERMAX

Maximum number of iterations, log_2(2/ACCURACY).

erradd

private double erradd

Constructor Detail

OrbitModel

public OrbitModel(VectorG[] times,
                  double[] y,
                  double[] err,
                  double off)

Measures at differnt point-in-time. Each VectorG has only a single dimentsion, this is the time.

Parameters:

times - The times of the measure, a Vector1D array.

y - The measurements of radial velocity.

err - The errors of the measurements.

Method Detail

getParameterCount

public int getParameterCount()

An orbit has six parameters, namely at indices:

• 0: The perihelium time, in the units of our time, like JD.
• 1: The orbital period, in our time units (e.g. days).
• 2: Excentricity of the orbit, dimensionless, 0-1.
• 3: The perihel argument, in radian.
• 4: The amplitude of the radial velocity, i.e. half of peak-to-valley, in units of the measurement.
• 5: The barycentric velocity of the system.

Returns:

The number of parameters in the model. This equals the dimension of the input vector in DataModel.getModel(vec_math.VectorG).

getMeasurementErrors

public double[] getMeasurementErrors()

If we had an external error, we subtract it again.

Specified by:

getMeasurementErrors in interface DataModel

Overrides:

getMeasurementErrors in class AbstractDataModel

getModel

public double[] getModel(VectorG a)

The model that fits an SB-1 orbit is given by

       v_r=v_0+e·K*cos(ω)+K*cos(ω+ν),
       

where v_0 is the baricentric velocity of the system, K is the radial velocity amplitude, e is the excentricity of the system, ω is the argument of perihelion, and ν is the true anomaly.

Specified by:

getModel in interface DataModel

Overrides:

getModel in class AbstractDataModel

Parameters:

a - The model parameters

Returns:

A model for the data points depending on the model parameters.

See Also:

getParameterCount

evaluateModel

public double evaluateModel(VectorG a,
                            VectorG t)

The model that fits an SB-1 orbit is given by

       v_r=v_0+e·K*cos(ω)+K*cos(ω+ν),
       

where v_0 is the baricentric velocity of the system, K is the radial velocity amplitude, e is the excentricity of the system, ω is the argument of perihelion, and ν is the true anomaly.

Parameters:

a - The fitted parameter

t - Only the zero'th component is used, this is the jd.

See Also:

getParameterCount

getOrbitSolution

private static double getOrbitSolution(double t,
                                       VectorG a)

Calculates the orbit rv for the given time and the specified parameter set.

Parameters:

t - The time, same units as in the model

a - The orbit model parameters, see getParameterCount() for order.

getGradientMatrix

public Matrix getGradientMatrix(VectorG a)

This is the derivative matrix of the orbital solution with respect to the parameters. From the definition of v_r one gets

Let v_nu be the partial derivative of vr with respect to ν

       v_nu = -K*sin(&omega+&nu);
       

Then, calling E the excentric anomaly, 2π/P*(t-t_0)=E-e·sinE, the respective gradients columns are:

• 0:
• 1:
• 2: Excentricity of the orbit, dimensionless, 0-1.
• 3:
• 4:
• 5: 1

getExcentricAnomaly

private static final double getExcentricAnomaly(double t,
                                                double t0,
                                                double p,
                                                double e)

From the time of perihelion passage, the orbital period, and the numeric eccentricity from the orbit, get the excentric anomaly for the current time. We solve Kepler's equation

       M=2π/P·(t-t0)=E-e·sin(E)
       

with Newton's method. If that does not converge, we can trap a solution between M+1 and M-1 and bisector from there.

getTrueAnomaly

private static final double getTrueAnomaly(double e,
                                           double exan)

From the eccentricity and the eccentric anomaly return the true anomaly.

       tan(ν/2) = sqrt((1+e)/(1-e))·tan(E/2)