List of usage examples for org.apache.commons.math3.util MathUtils TWO_PI
double TWO_PI
To view the source code for org.apache.commons.math3.util MathUtils TWO_PI.
Click Source Link
From source file:com.github.jonmarsh.waveform_processing_for_imagej.WaveformUtils.java
/** * Returns an array of window function values for a given length. Window * functions that require a parameter use the value {@code windowParam}. The * functions follow the asymmetric convention, in that the first value is * (typically) zero but the last value is nonzero. The window functions are * defined over the indices {@code i=0} to {@code i=n-1}, with the following * formulae:/* w ww . j a va2 s .c o m*/ * <ul> * <li> {@code BOHMAN: w[i]=((1-(abs(i-n/2)/(n/2)))*cos(PI*(abs(i-n/2)/(n/2))) + (1/PI)*sin(PI*(abs(i-n/2)/(n/2))), i=0,1,...,n-1 * } </li> * <li> {@code BLACKMAN: w[i]=0.42-0.5*cos(i*2*PI/n)+0.08*cos(2*i*2*PI/n), i=0,1,...,n-1 * } </li> * <li> {@code BLACKMAN_NUTTALL: w[i]=0.3635819-0.4891775*cos(i*2*PI/n)+0.1365995cos(2*i*2*PI/n)-0.0106411cos(3*i*2*PI/n), i=0,1,...,n-1 * } </li> * <li> {@code BLACKMAN_HARRIS: w[i]=0.42323-0.49755*cos(i*2*PI/n)+0.07922*cos(2*i*2*PI/n), i=0,1,...,n-1 * } </li> * <li> * {@code COSINE_TAPERED: w[i]=0.5*(1-cos(PI*i/m)), i=0,1,...,m-1; w[i]=0.5*(1-cos(PI*(n-i-1)/m)), i=n-m,...,n-1; w[i] = 1.0} * elsewhere; {@code m=floor(n*windowParam/2.0), 0.0<=windowParam<=1.0 } * </li> * <li> * {@code EXACT_BLACKMAN: w[i]=(7938-9240*cos(i*2*PI/n)+1430*cos(2*i*2*PI/n))/18608, i=0,1,...,n-1} * </li> * <li> * {@code EXPONENTIAL: w[i]=exp(a*i), a=ln(windowParam)/(n-1), i=0,1,...,n-1} * </li> * <li> * {@code FLAT_TOP: w[i]=0.21557895-0.41663158*cos(i*2*PI/n)+0.277263158*cos(2*i*2*PI/n)-0.083578947*cos(3*i*2*PI/n)+0.006947368*cos(4*i*2*PI/n), i=0,1,...,n-1} * </li> * <li> * {@code GAUSSIAN: w[i]=exp(-(i-(n/2))*(i-(n/2))/(2*windowParam*windowParam*(n+1)*(n+1))), i=0,1,...,n-1} * </li> * <li> {@code HAMMING: w[i]=0.54-0.46*cos(2*PI*i/n), i=0,1,...,n-1} </li> * <li> {@code HANNING: w[i]=0.5*(1.0-cos(2*PI*i/n)), i=0,1,...,n-1} </li> * <li> * {@code KAISER: w[i]=besselI0(windowParam*sqrt(1-a*a))/besselI0(windowParam), a=(i-k)/k, k=0.5*n, i=0,1,...,n-1, besselI0(x)} * is the modified Bessel function I<sub>0</sub>(x) </li> * <li> * {@code MODIFIED_BARTLETT_HANNING: w[i]=0.62-0.48*abs((i/n)-0.5)+0.38*cos(2*PI*((i/n)-0.5))}; * </li> * <li> * {@code PARZEN: w[i]=1.0-6.0*c*c+6.0*c*c*c for 0<=c<=0.5, and w[i]=2.0*(1.0-c)^3 for 0.5<c<1, c=abs(i-0.5*n)/(0.5*n), i=0,1,...,n-1} * </li> * <li> {@code TRIANGLE: w[i]=1.0-abs((2.0*i - n)/n), i=0,1,...,n-1} </li> * <li> * {@code WELCH: w[i] = 1.0-((i - 0.5*n)/(0.5*n))*((i - 0.5*n)/(0.5*n)), i=0,1,...,n-1} * </li> * </ul> * The default value is the unit array. * * @param windowType window function enumerated type constant as given * above * @param n length of window * @param windowParam used only for window functions that require it. * <ul> * <li> {@code COSINE_TAPERED}: {@code windowParam} refers to the ratio of * the length of the tapered section to the length of the entire signal, and * can have values from {@code 0.0} to {@code 1.0} </li> * <li> {@code EXPONENTIAL}: {@code windowParam} refers to the weighting * value of the last point of the window, and can be any nonzero positive * value </li> * <li> {@code GAUSSIAN}: {@code windowParam} refers to the * length-normalized standard deviation, and can be any number greater than * or equal to zero </li> * <li> {@code KAISER}: {@code windowParam} is proportional to the side-lobe * attenuation, and can be any real number </li> * </ul> * This value is ignored for other types of windows. * @param normalize set to true to normalize weights so that the sum of * all elements is equal to {@code 1.0} * @return an array of length {@code n} with weights computed via the * formulae above. If {@code n==1}, the returned array has a single * element equal to unity. If {@code n==0}, an empty array is * returned. {@code Null} is returned if {@code n<0}. If * {@code windowParam} is not within a valid range for a window * function that requires it, the return array is filled with zeros. */ public static final double[] windowFunction(WindowType windowType, int n, double windowParam, boolean normalize) { if (n > 1) { double[] w = new double[n]; switch (windowType) { case BLACKMAN: { double c = MathUtils.TWO_PI / n; for (int i = 0; i < n; i++) { w[i] = 0.42 - 0.5 * Math.cos(i * c) + 0.08 * Math.cos(2 * i * c); } break; } case BLACKMAN_HARRIS: { double c = MathUtils.TWO_PI / n; for (int i = 0; i < n; i++) { w[i] = 0.42323 - 0.49755 * Math.cos(i * c) + 0.07922 * Math.cos(2 * i * c); } break; } case BLACKMAN_NUTTALL: { double c = MathUtils.TWO_PI / n; for (int i = 0; i < n; i++) { w[i] = 0.3635819 - 0.4891775 * Math.cos(i * c) + 0.1365995 * Math.cos(2 * i * c) - 0.0106411 * Math.cos(3 * i * c); } break; } case BOHMAN: { double nOver2 = 0.5 * n; for (int i = 0; i < n; i++) { double fraction = Math.abs(i - nOver2) / nOver2; w[i] = ((1.0 - fraction) * Math.cos(Math.PI * fraction) + (1.0 / Math.PI) * (Math.sin(Math.PI * fraction))); } break; } case COSINE_TAPERED: { int m = (int) Math.floor(0.5 * n * windowParam); if (windowParam >= 0.0 && windowParam <= 1.0) { for (int i = 0; i < m; i++) { w[i] = 0.5 * (1.0 - Math.cos(Math.PI * i / (double) m)); } for (int i = m; i < n - m; i++) { w[i] = 1.0; } for (int i = n - m; i < n; i++) { w[i] = 0.5 * (1.0 - Math.cos(Math.PI * (n - i - 1) / (double) m)); } } break; } case EXACT_BLACKMAN: { double c = MathUtils.TWO_PI / n; for (int i = 0; i < n; i++) { w[i] = (7938.0 - 9240.0 * Math.cos(i * c) + 1430.0 * Math.cos(2 * i * c)) / 18608.0; } break; } case EXPONENTIAL: { if (windowParam > 0) { w[0] = 1.0; if (n > 1) { double a = Math.log(windowParam) / (n - 1.0); for (int i = 1; i < n; i++) { w[i] = Math.exp(a * i); } } } break; } case FLAT_TOP: { double a0 = 0.21557895, a1 = 0.41663158, a2 = 0.277263158, a3 = 0.083578947, a4 = 0.006947368; double c = MathUtils.TWO_PI / n; for (int i = 0; i < n; i++) { w[i] = a0 - a1 * Math.cos(i * c) + a2 * Math.cos(2 * i * c) - a3 * Math.cos(3 * i * c) + a4 * Math.cos(4 * i * c); } break; } case GAUSSIAN: { if (windowParam >= 0.0) { double m = 0.5 * n; double c = 1.0 / (2.0 * windowParam * windowParam * (n + 1) * (n + 1)); for (int i = 0; i < n; i++) { w[i] = Math.exp(-(i - m) * (i - m) * c); } } break; } case HAMMING: { double c = MathUtils.TWO_PI / n; for (int i = 0; i < n; i++) { w[i] = 0.54 - 0.46 * Math.cos(i * c); } break; } case HANNING: { double c = MathUtils.TWO_PI / n; for (int i = 0; i < n; i++) { w[i] = 0.5 * (1.0 - Math.cos(i * c)); } break; } case KAISER: { double k = 0.5 * n; for (int i = 0; i < n; i++) { double a = (i - k) / k; w[i] = besselI0(windowParam * Math.sqrt(1.0 - a * a)) / besselI0(windowParam); } break; } case MODIFIED_BARTLETT_HANNING: { for (int i = 0; i < n; i++) { double c = (double) i / (double) n - 0.5; w[i] = 0.62 - 0.48 * Math.abs(c) + 0.38 * Math.cos(MathUtils.TWO_PI * (c)); } break; } case PARZEN: { for (int i = 0; i < n; i++) { double c = Math.abs(i - 0.5 * n) / (0.5 * n); if (c <= 0.5) { w[i] = 1.0 - 6.0 * c * c + 6.0 * c * c * c; } else { w[i] = 2.0 * FastMath.pow(1.0 - c, 3); } } break; } case RECTANGLE: { for (int i = 0; i < n; i++) { w[i] = 1.0; } break; } case TRIANGLE: { for (int i = 0; i < n; i++) { w[i] = 1.0 - Math.abs((2.0 * i - n) / n); } break; } case WELCH: { for (int i = 0; i < n; i++) { double c = (i - 0.5 * n) / (0.5 * n); w[i] = 1.0 - c * c; } break; } default: { Arrays.fill(w, 1.0); } } if (normalize) { double sum = 0.0; for (int i = 0; i < n; i++) { sum += w[i]; } multiplyScalarInPlace(w, 1.0 / sum); } return w; } else if (n == 1) { return new double[] { 1.0 }; } else if (n == 0) { return new double[] {}; } else { return null; } }
From source file:com.github.jonmarsh.waveform_processing_for_imagej.WaveformUtils.java
/** * Returns an array of values for the desired window function of odd length * centered at index {@code 0}. The returned array is equivalent to positive * x-axis side of the symmetrical window function of length * {@code 2*radius+1} when the window is centered at the origin (i.e., when * the maximum value of {@code 1.0} is at index {@code 0} of the array). All * values in the returned array are nonzero; e.g., for a {@code TRIANGLE} * window function given an input of {@code radius=1}, the output is a * 2-element array the returned array is {@code {1.0, 0.5}} * (non-normalized)./*from w w w . java2 s . com*/ * * @param windowType window function enumerated type constant * @param radius number of nonzero-valued points for this window * function on the right side of the maximum value of * {@code 1.0} at index {@code 0}. * @param windowParam used only for window functions that require it. * <ul> * <li> {@code COSINE_TAPERED}: {@code windowParam} refers to the ratio of * the length of the tapered section to the length of the entire signal, and * can have values from {@code 0.0} to {@code 1.0} </li> * <li> {@code EXPONENTIAL}: {@code windowParam} refers to the weighting * value of the last point of the window, and can be any nonzero positive * value </li> * <li> {@code GAUSSIAN}: {@code windowParam} refers to the * length-normalized standard deviation, and can be any number greater than * or equal to zero </li> * <li> {@code KAISER}: {@code windowParam} is proportional to the side-lobe * attenuation, and can be any real number </li> * </ul> * This value is ignored for other types of windows. * @param normalize set to true to normalize weights so that the sum of * all elements of the equivalent <b>two-sided</b> window * function is equal to {@code 1.0} * @return an array of window function values of length {@code radius+1} * whose first element is always 1.0 and last element is nonzero. If * {@code radius==0}, the returned array has a single element with * value 1.0. If {@code radius<0}, the return value is {@code null}. */ public static final double[] windowFunctionSingleSided(WindowType windowType, int radius, double windowParam, boolean normalize) { if (radius > 0) { double[] w = new double[radius + 1]; w[0] = 1.0; int n = 2 * (radius + 1); switch (windowType) { case BLACKMAN: { double c = MathUtils.TWO_PI / n; for (int i = 1, j = i + radius + 1; i <= radius; i++, j++) { w[i] = 0.42 - 0.5 * Math.cos(j * c) + 0.08 * Math.cos(2 * j * c); } break; } case BLACKMAN_HARRIS: { double c = Math.PI / radius; for (int i = 1, j = i + radius; i <= radius; i++, j++) { w[i] = 0.42323 - 0.49755 * Math.cos(j * c) + 0.07922 * Math.cos(2 * j * c); } break; } case BLACKMAN_NUTTALL: { double c = Math.PI / radius; for (int i = 1, j = i + radius; i <= radius; i++, j++) { w[i] = 0.3635819 - 0.4891775 * Math.cos(j * c) + 0.1365995 * Math.cos(2 * j * c) - 0.0106411 * Math.cos(3 * j * c); } break; } case BOHMAN: { double nOver2 = 0.5 * n; for (int i = 1, j = i + radius + 1; i <= radius; i++, j++) { double fraction = Math.abs(j - nOver2) / nOver2; w[i] = ((1.0 - fraction) * Math.cos(Math.PI * fraction) + (1.0 / Math.PI) * (Math.sin(Math.PI * fraction))); } break; } case COSINE_TAPERED: { if (windowParam > 0.0 && windowParam <= 1.0) { int m = (int) Math.floor(radius * (1.0 - windowParam)); for (int i = 0; i < m; i++) { w[i] = 1.0; } for (int i = m; i <= radius; i++) { w[i] = 0.5 * (1.0 + Math.cos(Math.PI * (m - i) / (double) (radius + 1 - m))); } } break; } case EXACT_BLACKMAN: { double c = MathUtils.TWO_PI / (2 * radius); for (int i = 0, j = i + radius; i <= radius; i++, j++) { w[i] = (7938.0 - 9240.0 * Math.cos(j * c) + 1430.0 * Math.cos(2 * j * c)) / 18608.0; } break; } case EXPONENTIAL: { if (windowParam > 0) { double a = Math.log(windowParam) / radius; for (int i = 1; i <= radius; i++) { w[i] = Math.exp(a * i); } } break; } case FLAT_TOP: { double a0 = 0.21557895, a1 = 0.41663158, a2 = 0.277263158, a3 = 0.083578947, a4 = 0.006947368; double c = MathUtils.TWO_PI / (2 * radius); for (int i = 1, j = i + radius; i <= radius; i++, j++) { w[i] = a0 - a1 * Math.cos(j * c) + a2 * Math.cos(2 * j * c) - a3 * Math.cos(3 * j * c) + a4 * Math.cos(4 * j * c); } break; } case GAUSSIAN: { if (windowParam >= 0.0) { double c = 1.0 / (2.0 * windowParam * windowParam * (n - 1) * (n - 1)); for (int i = 1; i <= radius; i++) { w[i] = Math.exp(-(i * i) * c); } } break; } case HAMMING: { double c = Math.PI / radius; for (int i = 1, j = i + radius; i <= radius; i++, j++) { w[i] = 0.54 - 0.46 * Math.cos(j * c); } break; } case HANNING: { double c = Math.PI / (radius + 1); for (int i = 1, j = i + radius + 1; i <= radius; i++, j++) { w[i] = 0.5 * (1.0 - Math.cos(j * c)); } break; } case KAISER: { double k = radius; for (int i = 1, j = i + radius; i <= radius; i++, j++) { double a = (j - k) / k; w[i] = besselI0(windowParam * Math.sqrt(1.0 - a * a)) / besselI0(windowParam); } break; } case MODIFIED_BARTLETT_HANNING: { for (int i = 1, j = i + radius + 1; i <= radius; i++, j++) { double c = (double) j / (double) n - 0.5; w[i] = 0.62 - 0.48 * c + 0.38 * Math.cos(MathUtils.TWO_PI * (c)); } break; } case PARZEN: { for (int i = 1, j = i + radius + 1; i <= radius; i++, j++) { double c = Math.abs(j - 0.5 * n) / (0.5 * n); if (c <= 0.5) { w[i] = 1.0 - 6.0 * c * c + 6.0 * c * c * c; } else { w[i] = 2.0 * FastMath.pow(1.0 - c, 3); } } break; } case RECTANGLE: { for (int i = 1; i <= radius; i++) { w[i] = 1.0; } break; } case TRIANGLE: { for (int i = 1, j = i + radius + 1; i <= radius; i++, j++) { w[i] = 1.0 - Math.abs((2.0 * j - n) / n); } break; } case WELCH: { for (int i = 1, j = i + radius + 1; i <= radius; i++, j++) { double c = (j - 0.5 * n) / (0.5 * n); w[i] = 1.0 - c * c; } break; } default: { Arrays.fill(w, 1.0); } } if (normalize) { double sum = 1.0; for (int i = 1; i <= radius; i++) { sum += 2.0 * w[i]; } multiplyScalarInPlace(w, 1.0 / sum); } return w; } else if (radius == 0) { return new double[] { 1.0 }; } else { return null; } }
From source file:org.asoem.greyfish.core.space.WalledPointSpace.java
@Override public void moveObject(final O object, final Motion2D motion) { checkNotNull(object);/*from w w w .ja v a 2s . c o m*/ checkNotNull(motion); final double translation = motion.getTranslation(); final double rotation = motion.getRotation(); if (translation == 0 && rotation == 0) { return; } if (translation < 0) { throw new IllegalStateException("Translations < 0 are not supported: " + translation); } final double newOrientation = (rotation == 0) ? 0 : ((rotation) % MathUtils.TWO_PI + MathUtils.TWO_PI) % MathUtils.TWO_PI; final Point2D currentProjection = getProjection(object); final Point2D anchorPoint = checkNotNull(currentProjection, "Projection of {} is null", object) .getCentroid(); if (translation != 0) { final Point2D preferredPoint = ImmutablePoint2D.sum(anchorPoint, polarToCartesian(newOrientation, translation)); final Point2D maxPoint = maxTransition(anchorPoint, preferredPoint); //final MotionObject2D projection = MotionObject2DImpl.of(maxPoint.getX(), // maxPoint.getY(), !preferredPoint.equals(maxPoint)); //object.setProjection(projection); } // else: object.setProjection(MotionObject2DImpl.of(anchorPoint.getX(), // anchorPoint.getY(), newOrientation, false)); tree.invalidate(); }
From source file:org.asoem.greyfish.core.space.WalledPointSpace.java
@Override public double distance(final O agent, final double degrees) { checkNotNull(agent);//from ww w.j av a 2 s . c o m checkArgument(degrees >= 0 && degrees < MathUtils.TWO_PI, "Degrees must be in [0, TWO_PI), was %s", degrees); Point2D borderIntersection; final Point2D origin = getProjection(agent); if (origin == null) { throw new IllegalArgumentException("Has no projection: " + agent); } final ImmutablePoint2D destination = ImmutablePoint2D.sum(origin, polarToCartesian(degrees, Double.MAX_VALUE)); if (degrees < 90) { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), 0, 0, width(), 0); if (borderIntersection == null) { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), width(), 0, width(), height()); } } else if (degrees < 180) { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), width(), 0, width(), height()); if (borderIntersection == null) { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), 0, height(), width(), height()); } } else if (degrees < 270) { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), 0, height(), width(), height()); if (borderIntersection == null) { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), 0, 0, 0, height()); } } else { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), 0, 0, 0, height()); if (borderIntersection == null) { borderIntersection = intersection(origin.getX(), origin.getY(), destination.getX(), destination.getY(), 0, 0, width(), 0); } } assert borderIntersection != null; // There must always be an intersection with one border return origin.distance(maxTransition(origin, borderIntersection)); }
From source file:org.nd4j.linalg.api.rng.distribution.impl.SaddlePointExpansion.java
/** * Compute the logarithm of the PMF for a binomial distribution * using the saddle point expansion./*from w w w . j ava 2s .c om*/ * * @param x the value at which the probability is evaluated. * @param n the number of trials. * @param p the probability of success. * @param q the probability of failure (1 - p). * @return log(p(x)). */ public static double logBinomialProbability(int x, int n, double p, double q) { double ret; if (x == 0) { if (p < 0.1) { ret = -getDeviancePart(n, n * q) - n * p; } else { ret = n * FastMath.log(q); } } else if (x == n) { if (q < 0.1) { ret = -getDeviancePart(n, n * p) - n * q; } else { ret = n * FastMath.log(p); } } else { ret = getStirlingError(n) - getStirlingError(x) - getStirlingError(n - x) - getDeviancePart(x, n * p) - getDeviancePart(n - x, n * q); double f = (MathUtils.TWO_PI * x * (n - x)) / n; ret = -0.5 * FastMath.log(f) + ret; } return ret; }
From source file:org.orekit.forces.gravity.potential.ICGEMFormatReaderTest.java
private void checkValue(final double value, final AbsoluteDate date, final int n, final int m, final int refYear, final int refMonth, final int refDay, final double constant, final double trend, final double cosYear, final double sinYear, final double cosHalfYear, final double sinHalfYear, final int maxUlps) throws OrekitException { double factor = GravityFieldFactory.getUnnormalizationFactors(n, m)[n][m]; AbsoluteDate refDate = new AbsoluteDate(refYear, refMonth, refDay, 12, 0, 0, TimeScalesFactory.getTT()); double dtYear = date.durationFrom(refDate) / Constants.JULIAN_YEAR; double normalized = factor * (constant + trend * dtYear + cosYear * FastMath.cos(MathUtils.TWO_PI * dtYear) + sinYear * FastMath.sin(MathUtils.TWO_PI * dtYear) + cosHalfYear * FastMath.cos(MathUtils.TWO_PI * dtYear * 2) + sinHalfYear * FastMath.sin(MathUtils.TWO_PI * dtYear * 2)); double epsilon = maxUlps * FastMath.ulp(normalized); Assert.assertEquals(normalized, value, epsilon); }
From source file:org.orekit.forces.gravity.potential.PulsatingSphericalHarmonics.java
/** Simple constructor. * @param provider underlying part of the field * @param period period of the pulsation (s) * @param cosC cosine component of the cosine coefficients * @param sinC sine component of the cosine coefficients * @param cosS cosine component of the sine coefficients * @param sinS sine component of the sine coefficients *///from ww w .j a v a 2 s. co m PulsatingSphericalHarmonics(final RawSphericalHarmonicsProvider provider, final double period, final double[][] cosC, final double[][] sinC, final double[][] cosS, final double[][] sinS) { this.provider = provider; this.pulsation = MathUtils.TWO_PI / period; this.cosC = cosC; this.sinC = sinC; this.cosS = cosS; this.sinS = sinS; }
From source file:org.orekit.frames.TopocentricFrame.java
/** Get the azimuth of a point with regards to the topocentric frame center point. * <p>The azimuth is the angle between the North direction at local point and * the projection in local horizontal plane of the direction from local point * to given point. Azimuth angles are counted clockwise, i.e positive towards the East.</p> * @param extPoint point for which elevation shall be computed * @param frame frame in which the point is defined * @param date computation date//from ww w .j av a2 s. c o m * @return azimuth of the point * @exception OrekitException if frames transformations cannot be computed */ public double getAzimuth(final Vector3D extPoint, final Frame frame, final AbsoluteDate date) throws OrekitException { // Transform given point from given frame to topocentric frame final Transform t = getTransformTo(frame, date).getInverse(); final Vector3D extPointTopo = t.transformPosition(extPoint); // Compute azimuth double azimuth = FastMath.atan2(extPointTopo.getX(), extPointTopo.getY()); if (azimuth < 0.) { azimuth += MathUtils.TWO_PI; } return azimuth; }
From source file:org.orekit.frames.VEISProvider.java
/** Get the transform from GTOD at specified date. * @param date new value of the date// www . ja v a 2s . co m * @return transform at the specified date * @exception OrekitException if data embedded in the library cannot be read */ public Transform getTransform(final AbsoluteDate date) throws OrekitException { // offset from FIFTIES epoch (UT1 scale) final double dtai = date.durationFrom(VST_REFERENCE); final double dutc = TimeScalesFactory.getUTC().offsetFromTAI(date); final double dut1 = 0.0; // fixed at 0 since Veis parent is GTOD frame WITHOUT EOP corrections final double tut1 = dtai + dutc + dut1; final double ttd = tut1 / Constants.JULIAN_DAY; final double rdtt = ttd - (int) ttd; // compute Veis sidereal time, in radians final double vst = (VST0 + VST1 * ttd + MathUtils.TWO_PI * rdtt) % MathUtils.TWO_PI; // compute angular rotation of Earth, in rad/s final Vector3D rotationRate = new Vector3D(-VSTD, Vector3D.PLUS_K); // set up the transform from parent GTOD return new Transform(date, new Rotation(Vector3D.PLUS_K, vst), rotationRate); }
From source file:org.orekit.propagation.analytical.tle.SDP4.java
/** Computes SPACETRACK#3 compliant earth rotation angle. * @param date the current date//www . ja v a2 s . c o m * @return the ERA (rad) * @exception OrekitException when UTC time steps can't be read */ protected static double thetaG(final AbsoluteDate date) throws OrekitException { // Reference: The 1992 Astronomical Almanac, page B6. final double omega_E = 1.00273790934; final double jd = (date.durationFrom(AbsoluteDate.JULIAN_EPOCH) + date.timeScalesOffset(TimeScalesFactory.getUTC(), TimeScalesFactory.getTT())) / Constants.JULIAN_DAY; // Earth rotations per sidereal day (non-constant) final double UT = (jd + 0.5) % 1; final double seconds_per_day = Constants.JULIAN_DAY; final double jd_2000 = 2451545.0; /* 1.5 Jan 2000 = JD 2451545. */ final double t_cen = (jd - UT - jd_2000) / 36525.; double GMST = 24110.54841 + t_cen * (8640184.812866 + t_cen * (0.093104 - t_cen * 6.2E-6)); GMST = (GMST + seconds_per_day * omega_E * UT) % seconds_per_day; if (GMST < 0.) { GMST += seconds_per_day; } return MathUtils.TWO_PI * GMST / seconds_per_day; }