/* * HistMath.java * * Created on February 22, 2001, 1:12 PM */ package hep.aida.ref.histogram; import hep.aida.*; /** * * @author The AIDA team @ SLAC. * @version $Id: HistMath.java,v 1.9 2003/10/06 21:32:10 turri Exp $ */ class HistMath { static final double relPrec = 1E-10; protected HistMath() { } private double errorAdd(double e1, double e2) { return Math.sqrt(e1*e1+e2*e2); } private double errorSub(double e1, double e2) { return errorAdd(e1,e2); } private double errorMul(double e1, double l1, double e2, double l2) { return Math.sqrt(Math.pow(e1*l2, 2) + Math.pow(l1*e2, 2)); } private double errorDiv(double e1, double l1, double e2, double l2) { return Math.sqrt(Math.pow(e1/l2, 2) + Math.pow(e2*l1/(l2*l2), 2)); } private double addMean( double mean1, double height1, double mean2, double height2 ) { return ( mean1*height1 + mean2*height2 )/(height1+height2); } private double addRms( double rms1, double mean1, double height1, double rms2, double mean2, double height2 ) { double m = addMean( mean1,height1,mean2,height2); double h = height1+height2; double r = ((rms1*rms1*height1 + mean1*mean1*height1)+(rms2*rms2*height2 + mean2*mean2*height2))/h - m*m; if ( r < 0 ) { if ( Math.abs(r) < relPrec ) { r = 0; } else { r = 0; //throw new RuntimeException("Problem with rx "+r); } } return Math.sqrt(r); } private double subMean( double mean1, double height1, double mean2, double height2 ) { return ( mean1*height1 - mean2*height2 )/(height1-height2); } private double subRms( double rms1, double mean1, double height1, double rms2, double mean2, double height2 ) { double m = subMean( mean1,height1,mean2,height2); double h = height1-height2; double r = ((rms1*rms1*height1 + mean1*mean1*height1)-(rms2*rms2*height2 + mean2*mean2*height2))/h - m*m; if ( r < 0 ) { if ( Math.abs(r) < relPrec ) { r = 0; } else { r = 0; //throw new RuntimeException("Problem with r "+r); } } return Math.sqrt(r); } //----------------------------------------------------------------1D case------------------ /** *Checks for compatability of two axes */ private void checkCompatibility(IAxis axis1, IAxis axis2) throws IllegalArgumentException { if (!axis1.equals(axis2)) throw new IllegalArgumentException("Incompatible histogram binning"); } private void checkValidity(Histogram1D h1, Histogram1D h2) throws IllegalArgumentException { checkCompatibility(h1.axis(),h2.axis()); } /** * Adds two 1D Histogram * * @return h1 + h2 * @throws IllegalArgumentException if histogram binnings are incompatible */ IHistogram1D add(String name, IHistogram1D h1, IHistogram1D h2) throws IllegalArgumentException { return add( name, (Histogram1D)h1, (Histogram1D)h2 ); } IHistogram1D add(String name, Histogram1D h1, Histogram1D h2) throws IllegalArgumentException { checkValidity(h1, h2); int bins = h1.axis().bins()+2; double[] newHeights = new double[bins]; double[] newErrors = new double[bins]; double[] newMeans = new double[bins]; double[] newRmss = new double[bins]; int[] newEntries = new int [bins]; Histogram1D hist = new Histogram1D(name, name,h1.axis(),""); for(int i=IAxis.UNDERFLOW_BIN; i stop ) throw new IllegalArgumentException("Invalid indexes "+start+" "+stop); int xbins = h.xAxis().bins()+2; double[] newHeights = new double[xbins]; double[] newErrors = new double[xbins]; double[] newMeans = new double[xbins]; double[] newRmss = new double[xbins]; int[] newEntries = new int[xbins]; Histogram1D hist = new Histogram1D(name, name,h.xAxis()); for(int i=IAxis.UNDERFLOW_BIN; i stop ) throw new IllegalArgumentException("Invalid indexes "+start+" "+stop); int ybins = h.yAxis().bins()+2; double[] newHeights = new double[ybins]; double[] newErrors = new double[ybins]; double[] newMeans = new double[ybins]; double[] newRmss = new double[ybins]; int[] newEntries = new int[ybins]; Histogram1D hist = new Histogram1D(name, name,h.yAxis()); for(int j=IAxis.UNDERFLOW_BIN; j stop ) throw new IllegalArgumentException("Invalid indexes "+start+" "+stop); Histogram2D hist = new Histogram2D(name, name,h.xAxis(), h.yAxis()); int xbins = h.xAxis().bins()+2; int ybins = h.yAxis().bins()+2; double[][] newHeights = new double[xbins][ybins]; double[][] newErrors = new double[xbins][ybins]; double[][] newMeanXs = new double[xbins][ybins]; double[][] newRmsXs = new double[xbins][ybins]; double[][] newMeanYs = new double[xbins][ybins]; double[][] newRmsYs = new double[xbins][ybins]; int[][] newEntries = new int[xbins][ybins]; for(int i=IAxis.UNDERFLOW_BIN; i stop ) throw new IllegalArgumentException("Invalid indexes "+start+" "+stop); Histogram2D hist = new Histogram2D(name, name,h.yAxis(), h.zAxis()); int zbins = h.zAxis().bins()+2; int ybins = h.yAxis().bins()+2; double[][] newHeights = new double[ybins][zbins]; double[][] newErrors = new double[ybins][zbins]; double[][] newMeanYs = new double[ybins][zbins]; double[][] newRmsYs = new double[ybins][zbins]; double[][] newMeanZs = new double[ybins][zbins]; double[][] newRmsZs = new double[ybins][zbins]; int[][] newEntries = new int[ybins][zbins]; for(int j=IAxis.UNDERFLOW_BIN; j stop ) throw new IllegalArgumentException("Invalid indexes "+start+" "+stop); Histogram2D hist = new Histogram2D(name, name,h.xAxis(), h.zAxis()); int xbins = h.xAxis().bins()+2; int zbins = h.zAxis().bins()+2; double[][] newHeights = new double[xbins][zbins]; double[][] newErrors = new double[xbins][zbins]; double[][] newMeanXs = new double[xbins][zbins]; double[][] newRmsXs = new double[xbins][zbins]; double[][] newMeanZs = new double[xbins][zbins]; double[][] newRmsZs = new double[xbins][zbins]; int[][] newEntries = new int[xbins][zbins]; for(int i=IAxis.UNDERFLOW_BIN; i= 0 && index < nBins) return index+1; if (index == IAxis.UNDERFLOW_BIN) return 0; if (index == IAxis.OVERFLOW_BIN) return nBins+1; throw new IllegalArgumentException("Illegal argument "+index); } //Opposite conversion private int convertBackIndex( int index, int nBins ) { if ( index == 0 ) return IAxis.UNDERFLOW_BIN; if ( index == nBins+1 ) return IAxis.OVERFLOW_BIN; return index-1; } }