List of usage examples for java.lang Double isInfinite
public static boolean isInfinite(double v)
From source file:fr.amap.lidar.amapvox.commons.GTheta.java
/** * <p>Get the transmittance from the specified angle (radians or degrees) and the specified Leaf Angle Distribution.</p> * 0 is vertical, 90 is horizontal (zenithal angle, measured from vertical). * @param theta Angle//from ww w. ja v a2s .c om * @param degrees true if the given angle is in degrees, false otherwise * @return directional transmittance (GTheta) */ public double getGThetaFromAngle(double theta, boolean degrees) { if (degrees) { if (theta > 90) { //get an angle between 0 and 90 theta = 180 - theta; } } else { if (theta > (Math.PI / 2.0)) { //get an angle between 0 and pi/2 theta = Math.PI - theta; } } if (transmittanceFunctions != null && !isBuildingTable) { //a table was built int indice = 0; if (degrees) { indice = (int) (theta / res); } else { indice = (int) (Math.toDegrees(theta) / res); } if (indice >= transmittanceFunctions.length) { indice = transmittanceFunctions.length - 1; } else if (indice < 0) { indice = 0; } return transmittanceFunctions[indice]; } else { //no table was built, get transmittance on the fly if (pdfArray == null) { setupDensityProbabilityArray(DEFAULT_STEP_NUMBER); } if (distribution.getType() == SPHERIC) { return 0.5; //the result for spherical distribution is always 0.5, saving processing time } else { if (degrees) { theta = Math.toRadians(theta); } if (theta == 0) { theta = Double.MIN_VALUE; } if (theta >= Math.PI / 2.0) { theta = (Math.PI / 2.0) - 0.00001; } UnivariateFunction function1 = new CustomFunction1(theta); UnivariateFunction function2 = new CustomFunction2(theta); TrapezoidIntegrator integrator = new TrapezoidIntegrator(); double sum = 0; for (int j = 0; j < nbIntervals; j++) { double thetaL = (serie_angulaire[j] + serie_angulaire[j + 1]) / 2.0d; double Fi = (pdfArray[j]) / SOM; double cotcot = Math.abs(1 / (Math.tan(theta) * Math.tan(thetaL))); double Hi; if (cotcot > 1 || Double.isInfinite(cotcot)) { Hi = integrator.integrate(10000, function1, serie_angulaire[j], serie_angulaire[j + 1]); } else { Hi = integrator.integrate(10000, function2, serie_angulaire[j], serie_angulaire[j + 1]); //System.out.println("nb evaluations: " + integrator.getEvaluations()); } double Gi = Fi * Hi / ((Math.PI / 2) / (double) serie_angulaire.length); //because we need the average value not the actual integral value!!!! sum += Gi; } return sum; } } }
From source file:org.opentripplanner.api.resource.TimeGridWs.java
@GET @Produces({ "image/png" }) public Response getTimeGridPng(@QueryParam("base64") @DefaultValue("false") boolean base64) throws Exception { if (precisionMeters < 10) throw new IllegalArgumentException("Too small precisionMeters: " + precisionMeters); // Build the request RoutingRequest sptRequest = buildRequest(0); SampleGridRequest tgRequest = new SampleGridRequest(); tgRequest.maxTimeSec = maxTimeSec;//from w w w . ja v a 2 s .c om tgRequest.precisionMeters = precisionMeters; if (coordinateOrigin != null) tgRequest.coordinateOrigin = new GenericLocation(null, coordinateOrigin).getCoordinate(); // Get a sample grid ZSampleGrid<WTWD> sampleGrid = sampleGridRenderer.getSampleGrid(tgRequest, sptRequest); int cols = sampleGrid.getXMax() - sampleGrid.getXMin() + 1; int rows = sampleGrid.getYMax() - sampleGrid.getYMin() + 1; int channels = 4; // Hard-coded to RGBA // We force to 8 bits channel depth, some clients won't support more than 8 // (namely, HTML5 canvas...) ImageInfo imgInfo = new ImageInfo(cols, rows, 8, true, false, false); /** * TODO: PNGJ allow for progressive (ie line-by-line) writing. Thus we could theorically * prevent having to allocate a bit pixel array in the first place, but we would need a * line-by-line iterator on the sample grid, which is currently not the case. */ int[][] rgba = new int[rows][cols * channels]; ByteArrayOutputStream baos = new ByteArrayOutputStream(); PngWriter pw = new PngWriter(baos, imgInfo); pw.getMetadata().setText(PngChunkTextVar.KEY_Software, "OTPA"); pw.getMetadata().setText(PngChunkTextVar.KEY_Creation_Time, new Date().toString()); pw.getMetadata().setText(PngChunkTextVar.KEY_Description, "Sample grid bitmap"); String gridCornerStr = String.format(Locale.US, "%.8f,%.8f", sampleGrid.getCenter().y + sampleGrid.getYMin() * sampleGrid.getCellSize().y, sampleGrid.getCenter().x + sampleGrid.getXMin() * sampleGrid.getCellSize().x); String gridCellSzStr = String.format(Locale.US, "%.12f,%.12f", sampleGrid.getCellSize().y, sampleGrid.getCellSize().x); String offRoadDistStr = String.format(Locale.US, "%f", sampleGridRenderer.getOffRoadDistanceMeters(precisionMeters)); PngChunkTEXT gridCornerChunk = new PngChunkTEXT(imgInfo); gridCornerChunk.setKeyVal(OTPA_GRID_CORNER, gridCornerStr); pw.getChunksList().queue(gridCornerChunk); PngChunkTEXT gridCellSzChunk = new PngChunkTEXT(imgInfo); gridCellSzChunk.setKeyVal(OTPA_GRID_CELL_SIZE, gridCellSzStr); pw.getChunksList().queue(gridCellSzChunk); PngChunkTEXT offRoadDistChunk = new PngChunkTEXT(imgInfo); offRoadDistChunk.setKeyVal(OTPA_OFFROAD_DIST, offRoadDistStr); pw.getChunksList().queue(offRoadDistChunk); double unit; switch (zDataType) { case TIME: unit = 1.0; // 1:1sec, max 18h break; case BOARDINGS: unit = 1000.0; // 1:0.001 boarding, max 65.5 break; case WALK_DISTANCE: unit = 10.0; // 1:0.1m, max 6.55km break; default: throw new IllegalArgumentException("Unsupported Z DataType."); } for (ZSamplePoint<WTWD> p : sampleGrid) { WTWD z = p.getZ(); int row = p.getY() - sampleGrid.getYMin(); int col = p.getX() - sampleGrid.getXMin(); double zz; switch (zDataType) { case TIME: zz = z.wTime / z.w; break; case BOARDINGS: zz = z.wBoardings / z.w; break; case WALK_DISTANCE: zz = z.wWalkDist / z.w; break; default: throw new IllegalArgumentException("Unsupported Z DataType."); } int iz; if (Double.isInfinite(zz)) { iz = 65535; } else { iz = ImageLineHelper.clampTo_0_65535((int) Math.round(zz * unit)); if (iz == 65535) iz = 65534; // Clamp } // d is expressed as a percentage of grid size, max 255%. // Sometimes d will be bigger than 2.55 x grid size, // but this should not be too much important as we are off-bounds. int id = ImageLineHelper.clampTo_0_255((int) Math.round(z.d / precisionMeters * 100)); int offset = col * channels; rgba[row][offset + 0] = (iz & 0xFF); // z low 8 bits rgba[row][offset + 1] = (iz >> 8); // z high 8 bits rgba[row][offset + 2] = id; // d /* * Keep the alpha channel at 255, otherwise the RGB channel will be downsampled on some * rendering clients (namely, JS canvas). */ rgba[row][offset + 3] = 255; } for (int row = 0; row < rgba.length; row++) { ImageLineInt iline = new ImageLineInt(imgInfo, rgba[row]); pw.writeRow(iline, row); } pw.end(); // Disallow caching on client side CacheControl cc = new CacheControl(); cc.setNoCache(true); // Also put the meta-data in the HTML header (easier to read from JS) byte[] data = baos.toByteArray(); if (base64) { data = Base64.encodeBase64(data); } return Response.ok().cacheControl(cc).entity(data).header(OTPA_GRID_CORNER, gridCornerStr) .header(OTPA_GRID_CELL_SIZE, gridCellSzStr).header(OTPA_OFFROAD_DIST, offRoadDistStr).build(); }
From source file:com.opengamma.analytics.math.function.PiecewisePolynomialWithSensitivityFunction1D.java
/** * @param pp {@link PiecewisePolynomialResultsWithSensitivity} * @param xKey //from w w w .ja v a 2s .c om * @return Node sensitivity of second derivative value at x=xKey */ public DoubleMatrix1D differentiateTwiceNodeSensitivity(final PiecewisePolynomialResultsWithSensitivity pp, final double xKey) { ArgumentChecker.notNull(pp, "null pp"); ArgumentChecker.isFalse(Double.isNaN(xKey), "xKey containing NaN"); ArgumentChecker.isFalse(Double.isInfinite(xKey), "xKey containing Infinity"); if (pp.getDimensions() > 1) { throw new NotImplementedException(); } final int nCoefs = pp.getOrder(); ArgumentChecker.isFalse(nCoefs < 3, "Polynomial degree is too low"); final double[] knots = pp.getKnots().getData(); final int nKnots = knots.length; int interval = FunctionUtils.getLowerBoundIndex(knots, xKey); if (interval == nKnots - 1) { interval--; // there is 1 less interval that knots } final double s = xKey - knots[interval]; final DoubleMatrix2D a = pp.getCoefficientSensitivity(interval); DoubleMatrix1D res = (DoubleMatrix1D) MA.scale(a.getRowVector(0), (nCoefs - 1) * (nCoefs - 2)); for (int i = 1; i < nCoefs - 2; i++) { res = (DoubleMatrix1D) MA.scale(res, s); res = (DoubleMatrix1D) MA.add(res, MA.scale(a.getRowVector(i), (nCoefs - 1 - i) * (nCoefs - 2 - i))); } return res; }
From source file:com.opengamma.analytics.math.interpolation.MonotoneConvexSplineInterpolator.java
@Override public DoubleMatrix1D interpolate(final double[] xValues, final double[] yValues, final double[] x) { ArgumentChecker.notNull(x, "x"); final int keyLength = x.length; double[] res = new double[keyLength]; final PiecewisePolynomialResult result = interpolate(xValues, yValues); final DoubleMatrix2D coefsMatrixIntegrate = result.getCoefMatrix(); final int nKnots = coefsMatrixIntegrate.getNumberOfRows() + 1; final double[] knots = result.getKnots().getData(); for (int j = 0; j < keyLength; ++j) { int indicator = 0; if (x[j] <= knots[1]) { indicator = 0;/*from w w w . j a v a 2 s. c o m*/ } else { for (int i = 1; i < nKnots - 1; ++i) { if (knots[i] < x[j]) { indicator = i; } } } final double[] coefs = coefsMatrixIntegrate.getRowVector(indicator).getData(); res[j] = getValue(coefs, x[j], knots[indicator]); ArgumentChecker.isFalse(Double.isInfinite(res[j]), "Too large/small data values or xKey"); ArgumentChecker.isFalse(Double.isNaN(res[j]), "Too large/small data values or xKey"); } return new DoubleMatrix1D(res); }
From source file:ec.ui.view.res.ResidualsView.java
private Range calcRange(double[] values) { double min = Double.NEGATIVE_INFINITY, max = -Double.POSITIVE_INFINITY; DescriptiveStatistics stats = new DescriptiveStatistics(new DataBlock(values)); double smin = stats.getMin(), smax = stats.getMax(); if (Double.isInfinite(min) || smin < min) { min = smin;/* ww w. j a v a 2s. co m*/ } if (Double.isInfinite(max) || smax > max) { max = smax; } if (Double.isInfinite(max) || Double.isInfinite(min)) { return new Range(0, 1); } double length = max - min; if (length == 0) { return new Range(0, 1); } else { //double eps = length * .05; //return new Range(min - eps, max + eps); return new Range(min, max); } }
From source file:org.nd4j.linalg.solvers.VectorizedBackTrackLineSearch.java
public double optimize(INDArray line, int lineSearchIteration, double initialStep) throws InvalidStepException { INDArray g, x, oldParameters;//from ww w. ja v a2 s . c om double slope, test, alamin, alam, alam2, tmplam; double rhs1, rhs2, a, b, disc, oldAlam; double f, fold, f2; g = function.getValueGradient(lineSearchIteration); // gradient x = function.getParameters(); // parameters oldParameters = x.dup(); alam2 = tmplam = 0.0f; f2 = fold = function.getValue(); if (logger.isDebugEnabled()) { logger.trace("ENTERING BACKTRACK\n"); logger.trace("Entering BackTrackLnSrch, value=" + fold + ",\ndirection.oneNorm:" + line.norm1(Integer.MAX_VALUE) + " direction.infNorm:" + FastMath.max(Float.NEGATIVE_INFINITY, (double) Transforms.abs(line).max(Integer.MAX_VALUE).element())); } LinAlgExceptions.assertValidNum(g); double sum = (double) line.norm2(Integer.MAX_VALUE).element(); if (sum > stpmax) { logger.warn("attempted step too big. scaling: sum= " + sum + ", stpmax= " + stpmax); line.muli(stpmax / sum); } //dot product slope = Nd4j.getBlasWrapper().dot(g, line); logger.debug("slope = " + slope); if (slope < 0) { throw new InvalidStepException("Slope = " + slope + " is negative"); } if (slope == 0) throw new InvalidStepException("Slope = " + slope + " is zero"); // find maximum lambda // converge when (delta x) / x < REL_TOLX for all coordinates. // the largest step size that triggers this threshold is // precomputed and saved in alamin INDArray maxOldParams = Nd4j.create(line.length()); for (int i = 0; i < line.length(); i++) { maxOldParams.putScalar(i, Math.max(Math.abs(oldParameters.getDouble(i)), 1.0)); } INDArray testMatrix = Transforms.abs(line).div(maxOldParams); test = testMatrix.max(Integer.MAX_VALUE).getDouble(0); alamin = relTolx / test; alam = 1.0f; oldAlam = 0.0f; int iteration = 0; // look for step size in direction given by "line" for (iteration = 0; iteration < maxIterations; iteration++) { function.setCurrentIteration(lineSearchIteration); // x = oldParameters + alam*line // initially, alam = 1.0, i.e. take full Newton step logger.trace("BackTrack loop iteration " + iteration + " : alam=" + alam + " oldAlam=" + oldAlam); logger.trace("before step, x.1norm: " + x.norm1(Integer.MAX_VALUE) + "\nalam: " + alam + "\noldAlam: " + oldAlam); assert (alam != oldAlam) : "alam == oldAlam"; x.addi(line.mul(alam - oldAlam)); // step double norm1 = x.norm1(Integer.MAX_VALUE).getDouble(0); logger.debug("after step, x.1norm: " + norm1); // check for convergence //convergence on delta x if ((alam < alamin) || smallAbsDiff(oldParameters, x)) { // if ((alam < alamin)) { function.setParameters(oldParameters); f = function.getValue(); logger.trace("EXITING BACKTRACK: Jump too small (alamin=" + alamin + "). Exiting and using xold. Value=" + f); return 0.0f; } function.setParameters(x); oldAlam = alam; f = function.getValue(); logger.debug("value = " + f); // sufficient function increase (Wolf condition) if (f >= fold + ALF * alam * slope) { logger.debug("EXITING BACKTRACK: value=" + f); if (f < fold) throw new IllegalStateException("Function did not increase: f=" + f + " < " + fold + "=fold"); return alam; } // if value is infinite, i.e. we've // jumped to unstable territory, then scale down jump else if (Double.isInfinite(f) || Double.isInfinite(f2)) { logger.warn("Value is infinite after jump " + oldAlam + ". f=" + f + ", f2=" + f2 + ". Scaling back step size..."); tmplam = .2f * alam; if (alam < alamin) { //convergence on delta x function.setParameters(oldParameters); f = function.getValue(); logger.warn("EXITING BACKTRACK: Jump too small. Exiting and using xold. Value=" + f); return 0.0f; } } else { // backtrack if (alam == 1.0) // first time through tmplam = -slope / (2.0f * (f - fold - slope)); else { rhs1 = f - fold - alam * slope; rhs2 = f2 - fold - alam2 * slope; assert ((alam - alam2) != 0) : "FAILURE: dividing by alam-alam2. alam=" + alam; a = (rhs1 / (double) (FastMath.pow(alam, 2)) - rhs2 / (double) (FastMath.pow(alam2, 2))) / (alam - alam2); b = (-alam2 * rhs1 / (alam * alam) + alam * rhs2 / (alam2 * alam2)) / (alam - alam2); if (a == 0.0) tmplam = -slope / (2.0f * b); else { disc = b * b - 3.0f * a * slope; if (disc < 0.0) { tmplam = .5f * alam; } else if (b <= 0.0) tmplam = (-b + (double) FastMath.sqrt(disc)) / (3.0f * a); else tmplam = -slope / (b + (double) FastMath.sqrt(disc)); } if (tmplam > .5f * alam) tmplam = .5f * alam; // lambda <= .5 lambda_1 } } alam2 = alam; f2 = f; logger.debug("tmplam:" + tmplam); alam = Math.max(tmplam, .1f * alam); // lambda >= .1*Lambda_1 } if (iteration >= maxIterations) throw new IllegalStateException("Too many iterations."); return 0.0f; }
From source file:org.deeplearning4j.optimize.solvers.VectorizedBackTrackLineSearch.java
public double optimize(INDArray line, int lineSearchIteration, double initialStep) throws InvalidStepException { INDArray g, x, oldParameters;//from ww w. j a v a2 s . com double slope, test, alamin, alam, alam2, tmplam; double rhs1, rhs2, a, b, disc, oldAlam; double f, fold, f2; g = function.getValueGradient(lineSearchIteration); // gradient x = function.getParameters(); // parameters oldParameters = x.dup(); alam2 = tmplam = 0.0f; f2 = fold = function.getValue(); if (logger.isDebugEnabled()) { logger.trace("ENTERING BACKTRACK\n"); logger.trace("Entering BackTrackLnSrch, value=" + fold + ",\ndirection.oneNorm:" + line.norm1(Integer.MAX_VALUE) + " direction.infNorm:" + FastMath.max(Float.NEGATIVE_INFINITY, Transforms.abs(line).max(Integer.MAX_VALUE).getDouble(0))); } LinAlgExceptions.assertValidNum(g); double sum = line.norm2(Integer.MAX_VALUE).getDouble(0); if (sum > stpmax) { logger.warn("attempted step too big. scaling: sum= " + sum + ", stpmax= " + stpmax); line.muli(stpmax / sum); } //dot product slope = Nd4j.getBlasWrapper().dot(g, line); logger.debug("slope = " + slope); if (slope < 0) { throw new InvalidStepException("Slope = " + slope + " is negative"); } if (slope == 0) throw new InvalidStepException("Slope = " + slope + " is zero"); // find maximum lambda // converge when (delta x) / x < REL_TOLX for all coordinates. // the largest step size that triggers this threshold is // precomputed and saved in alamin INDArray maxOldParams = Nd4j.create(line.length()); for (int i = 0; i < line.length(); i++) { maxOldParams.putScalar(i, Math.max(Math.abs(oldParameters.getDouble(i)), 1.0)); } INDArray testMatrix = Transforms.abs(line).div(maxOldParams); test = testMatrix.max(Integer.MAX_VALUE).getDouble(0); alamin = relTolx / test; alam = 1.0f; oldAlam = 0.0f; int iteration = 0; // look for step size in direction given by "line" for (iteration = 0; iteration < maxIterations; iteration++) { function.setCurrentIteration(lineSearchIteration); // x = oldParameters + alam*line // initially, alam = 1.0, i.e. take full Newton step logger.trace("BackTrack loop iteration " + iteration + " : alam=" + alam + " oldAlam=" + oldAlam); logger.trace("before step, x.1norm: " + x.norm1(Integer.MAX_VALUE) + "\nalam: " + alam + "\noldAlam: " + oldAlam); assert (alam != oldAlam) : "alam == oldAlam"; x.addi(line.mul(alam - oldAlam)); // step double norm1 = x.norm1(Integer.MAX_VALUE).getDouble(0); logger.debug("after step, x.1norm: " + norm1); // check for convergence //convergence on delta x if ((alam < alamin) || smallAbsDiff(oldParameters, x)) { // if ((alam < alamin)) { function.setParameters(oldParameters); f = function.getValue(); logger.trace("EXITING BACKTRACK: Jump too small (alamin=" + alamin + "). Exiting and using xold. Value=" + f); return 0.0f; } function.setParameters(x); oldAlam = alam; f = function.getValue(); logger.debug("value = " + f); // sufficient function increase (Wolf condition) if (f >= fold + ALF * alam * slope) { logger.debug("EXITING BACKTRACK: value=" + f); if (f < fold) throw new IllegalStateException("Function did not increase: f=" + f + " < " + fold + "=fold"); return alam; } // if value is infinite, i.e. we've // jumped to unstable territory, then scale down jump else if (Double.isInfinite(f) || Double.isInfinite(f2)) { logger.warn("Value is infinite after jump " + oldAlam + ". f=" + f + ", f2=" + f2 + ". Scaling back step size..."); tmplam = .2f * alam; if (alam < alamin) { //convergence on delta x function.setParameters(oldParameters); f = function.getValue(); logger.warn("EXITING BACKTRACK: Jump too small. Exiting and using xold. Value=" + f); return 0.0f; } } else { // backtrack if (alam == 1.0) // first time through tmplam = -slope / (2.0f * (f - fold - slope)); else { rhs1 = f - fold - alam * slope; rhs2 = f2 - fold - alam2 * slope; assert ((alam - alam2) != 0) : "FAILURE: dividing by alam-alam2. alam=" + alam; a = (rhs1 / (double) (FastMath.pow(alam, 2)) - rhs2 / (double) (FastMath.pow(alam2, 2))) / (alam - alam2); b = (-alam2 * rhs1 / (alam * alam) + alam * rhs2 / (alam2 * alam2)) / (alam - alam2); if (a == 0.0) tmplam = -slope / (2.0f * b); else { disc = b * b - 3.0f * a * slope; if (disc < 0.0) { tmplam = .5f * alam; } else if (b <= 0.0) tmplam = (-b + (double) FastMath.sqrt(disc)) / (3.0f * a); else tmplam = -slope / (b + (double) FastMath.sqrt(disc)); } if (tmplam > .5f * alam) tmplam = .5f * alam; // lambda <= .5 lambda_1 } } alam2 = alam; f2 = f; logger.debug("tmplam:" + tmplam); alam = Math.max(tmplam, .1f * alam); // lambda >= .1*Lambda_1 } if (iteration >= maxIterations) throw new IllegalStateException("Too many iterations."); return 0.0f; }
From source file:bide.math.NormalDistribution.java
/** A more accurate and faster implementation of the cdf (taken from function pnorm in the R statistical language) * This implementation has discrepancies depending on the programming language and system architecture * In Java, returned values become zero once z reaches -37.5193 exactly on the machine tested * In the other implementation, the returned value 0 at about z = -8 * In C, this 0 value is reached approximately z = -37.51938 * * Will later need to be optimised for BEAST * * @param x argument/* w w w . ja v a2 s. c o m*/ * @param mu mean * @param sigma standard deviation * @param log_p is p logged * @return cdf at x */ public static double cdf(double x, double mu, double sigma, boolean log_p) { boolean i_tail = false; double p, cp = Double.NaN; if (Double.isNaN(x) || Double.isNaN(mu) || Double.isNaN(sigma)) { return Double.NaN; } if (Double.isInfinite(x) && mu == x) { /* x-mu is NaN */ return Double.NaN; } if (sigma <= 0) { if (sigma < 0) { return Double.NaN; } return (x < mu) ? 0.0 : 1.0; } p = (x - mu) / sigma; if (Double.isInfinite(p)) { return (x < mu) ? 0.0 : 1.0; } x = p; if (Double.isNaN(x)) { return Double.NaN; } double xden, xnum, temp, del, eps, xsq, y; int i; boolean lower, upper; eps = DBL_EPSILON * 0.5; lower = !i_tail; upper = i_tail; y = Math.abs(x); if (y <= 0.67448975) { /* Normal.quantile(3/4, 1, 0) = 0.67448975 */ if (y > eps) { xsq = x * x; xnum = a[4] * xsq; xden = xsq; for (i = 0; i < 3; i++) { xnum = (xnum + a[i]) * xsq; xden = (xden + b[i]) * xsq; } } else { xnum = xden = 0.0; } temp = x * (xnum + a[3]) / (xden + b[3]); if (lower) { p = 0.5 + temp; } if (upper) { cp = 0.5 - temp; } if (log_p) { if (lower) { p = Math.log(p); } if (upper) { cp = Math.log(cp); } } } else if (y <= M_SQRT_32) { /* Evaluate pnorm for 0.67448975 = Normal.quantile(3/4, 1, 0) < |x| <= sqrt(32) ~= 5.657 */ xnum = c[8] * y; xden = y; for (i = 0; i < 7; i++) { xnum = (xnum + c[i]) * y; xden = (xden + d[i]) * y; } temp = (xnum + c[7]) / (xden + d[7]); //do_del(y); //swap_tail; //#define do_del(X) \ xsq = ((int) (y * CUTOFF)) * 1.0 / CUTOFF; del = (y - xsq) * (y + xsq); if (log_p) { p = (-xsq * xsq * 0.5) + (-del * 0.5) + Math.log(temp); if ((lower && x > 0.0) || (upper && x <= 0.0)) { cp = Math.log(1.0 - Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp); } } else { p = Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp; cp = 1.0 - p; } //#define swap_tail \ if (x > 0.0) { temp = p; if (lower) { p = cp; } cp = temp; } } /* else |x| > sqrt(32) = 5.657 : * the next two case differentiations were really for lower=T, log=F * Particularly *not* for log_p ! * Cody had (-37.5193 < x && x < 8.2924) ; R originally had y < 50 * Note that we do want symmetry(0), lower/upper -> hence use y */ else if (log_p || (lower && -37.5193 < x && x < 8.2924) || (upper && -8.2924 < x && x < 37.5193)) { /* Evaluate pnorm for x in (-37.5, -5.657) union (5.657, 37.5) */ xsq = 1.0 / (x * x); xnum = p_[5] * xsq; xden = xsq; for (i = 0; i < 4; i++) { xnum = (xnum + p_[i]) * xsq; xden = (xden + q[i]) * xsq; } temp = xsq * (xnum + p_[4]) / (xden + q[4]); temp = (M_1_SQRT_2PI - temp) / y; //do_del(x); xsq = ((int) (x * CUTOFF)) * 1.0 / CUTOFF; del = (x - xsq) * (x + xsq); if (log_p) { p = (-xsq * xsq * 0.5) + (-del * 0.5) + Math.log(temp); if ((lower && x > 0.0) || (upper && x <= 0.0)) { cp = Math.log(1.0 - Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp); } } else { p = Math.exp(-xsq * xsq * 0.5) * Math.exp(-del * 0.5) * temp; cp = 1.0 - p; } //swap_tail; if (x > 0.0) { temp = p; if (lower) { p = cp; } cp = temp; } } else { /* no log_p , large x such that probs are 0 or 1 */ if (x > 0) { p = 1.0; cp = 0.0; } else { p = 0.0; cp = 1.0; } } return p; }
From source file:org.structnetalign.ReportGenerator.java
/** * Yes, I really am that OCD.//from w ww . ja v a2s . c o m */ private static Object format(Object value) { if (value == null) return null; if (value.getClass().isAssignableFrom(Double.class)) { double x = Math.abs((double) value); String minus = ""; if (x < 0) minus = ""; if (Double.isInfinite(x)) { value = minus + ""; } else { value = minus + x; } } return value; }
From source file:org.eclipse.january.metadata.internal.StatisticsMetadataImpl.java
/** * Calculate summary statistics for a dataset * @param mm//www. j av a 2 s. c o m * @param ignoreNaNs if true, ignore NaNs * @param ignoreInfs if true, ignore infinities */ @SuppressWarnings("unchecked") private void setMaxMin(final MaxMin<T> mm, final boolean ignoreNaNs, final boolean ignoreInfs) { final IndexIterator iter = dataset.getIterator(); if (!DTypeUtils.isDTypeNumerical(dtype)) { // TODO FIXME for strings // treat non-numerical as strings in lexicographic order String smax = dataset.getStringAbs(0); String smin = smax; while (iter.hasNext()) { final String val = dataset.getStringAbs(iter.index); hash = (int) (hash * 19 + val.hashCode()); if (val.compareTo(smax) > 0) { smax = val; } if (val.compareTo(smin) < 0) { smin = val; } } hash = hash * 19 + dtype * 17 + isize; mm.maximum = (T) smax; mm.minimum = (T) smin; mm.maximumPositions = null; mm.minimumPositions = null; return; } if (isize == 1) { double amax = Double.NEGATIVE_INFINITY; double amin = Double.POSITIVE_INFINITY; boolean hasNaNs = false; if (dataset.hasFloatingPointElements() && (ignoreNaNs || ignoreInfs)) { while (iter.hasNext()) { final double val = dataset.getElementDoubleAbs(iter.index); hash = (int) (hash * 19 + Double.doubleToRawLongBits(val)); if (Double.isNaN(val)) { if (ignoreNaNs) continue; hasNaNs = true; } else if (Double.isInfinite(val)) { if (ignoreInfs) continue; } if (val > amax) { amax = val; } if (val < amin) { amin = val; } } } else if (dataset.hasFloatingPointElements()) { while (iter.hasNext()) { final double val = dataset.getElementDoubleAbs(iter.index); hash = (int) (hash * 19 + Double.doubleToRawLongBits(val)); if (Double.isNaN(val)) { hasNaNs = true; continue; } if (val > amax) { amax = val; } if (val < amin) { amin = val; } } } else { while (iter.hasNext()) { final long val = dataset.getElementLongAbs(iter.index); hash = (int) (hash * 19 + val); if (val > amax) { amax = val; } if (val < amin) { amin = val; } } } mm.maximum = (T) (hasNaNs ? Double.NaN : DTypeUtils.fromDoubleToBiggestNumber(amax, dtype)); mm.minimum = (T) (hasNaNs ? Double.NaN : DTypeUtils.fromDoubleToBiggestNumber(amin, dtype)); } else { while (iter.hasNext()) { for (int j = 0; j < isize; j++) { final double val = dataset.getElementDoubleAbs(iter.index + j); hash = (int) (hash * 19 + Double.doubleToRawLongBits(val)); } } mm.maximum = null; mm.minimum = null; } hash = hash * 19 + dtype * 17 + isize; mm.maximumPositions = null; mm.minimumPositions = null; }