newdynamicdataanalyzer.java

一个用于排队系统仿真的开源软件,有非常形象的图象仿真过程!

/**    
  * Copyright (C) 2006, Laboratorio di Valutazione delle Prestazioni - Politecnico di Milano

  * This program is free software; you can redistribute it and/or modify
  * it under the terms of the GNU General Public License as published by
  * the Free Software Foundation; either version 2 of the License, or
  * (at your option) any later version.

  * This program is distributed in the hope that it will be useful,
  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  * GNU General Public License for more details.

  * You should have received a copy of the GNU General Public License
  * along with this program; if not, write to the Free Software
  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
  */
  
package jmt.engine.dataAnalysis;

import Jama.Matrix;
import jmt.engine.math.BatchCircularList;
import jmt.engine.math.Spectrum;
import jmt.engine.math.TStudent;
import jmt.engine.math.WeightedMeanVar;

/**
 * <p>Title: New Dynamic Data Analyzer Class</p>
 * <p>Description: This class provides a new dynamic data analyzer method based
 * on Fisherman Rule, MSER-m Rule and Welch's Spectral test.</p>
 *
 * <p>Old one was slow, wrongly implemented and worked by "black magic".</p>
 *
 * @author Bertoli Marco (Fisman rule),
 *         Omini Stefano (null test),
 *         Granata Federico (spectral test)
 *
 *         Date: 11-gen-2006
 *         Time: 12.02.04
 */
public class NewDynamicDataAnalyzer implements DynamicDataAnalyzer{
    private double alpha, precision;
    private int maxData;

    // end tells if computation is ended
    private boolean end;

    // success tells if measure was calculated correctly or failed
    private boolean success = false;

    // Number of analyzed samples
    private int nsamples = 0;

    // Tells if transient has been deleted
    private boolean delet = false;

// ---- Variables for Fishman Means Rule -------------------------------------------------------
    // Number of crosses to be checked by Fishman Rule
    private static final int crossesNum = 19;
    // Number of samples to be collected to check Fishman Rule
    private static final int checkSize = 500;
    // Multiplier to space to be added when arrays are full
    private static final int hArrayMult = 4;
    // An epsilon added to the mean value to avoid problems with monothonic
    // measures tending to an asympthotic value
    private static final double hEpsilon = 1e-8;

    // Tells if heuristic was passed
    private boolean heuristicPassed = false;

    // Number of analyzed samples within this rule
    private int hSamplesNum = 0;

    // Arrays with mean values
    private double[] hMeans;

    // Transient Length
    private int transientLen = 0;
// ---------------------------------------------------------------------------------------------

// ---- Variables for MSER-m Rule --------------------------------------------------------------
    // Number of samples in a batch. 5 is raccomanded as will perform a MSER-5 rule that seems
    // the most stable
    private static final int batchSize = 5;

    // Maximum number of samples to be checked by the rule.
    // Big values are better but will REALLY slow simulation down
    private static final int maxRuleSamples = 5000;

    // Used to store transient datas to perform MSER-m rule and initialize spectral analisys
    private BatchCircularList batches = new BatchCircularList(batchSize, maxRuleSamples);

    /**minimum length (number of samples) of transient period */
    int minSamples = 0;
// ---------------------------------------------------------------------------------------------


// ---- Spectral Analysis ----------------------------------------------------------------------
    // Number of batches
    private int numBatch = 128;
    // Number of samples per batch
    private int batchLen = 8;

    // Mean of batches
    private double[] batchMean;
    // Mean of batches weighted
    private double[] weightBatchMean;

    /**index of the batch*/
    int batch;

    /**the extimated mean */
    protected double extMean = 0;
    /**the extimated variance */
    protected double extVar = 0;
    /**the exitmated confidence intervals */
    protected double confInt = 0;

    /**constant chosen to reduce bias of variance extimator */
    protected double C1;
    /**degrees of freeedom of the t-distribution */
    protected int C2;

    /**number of points of the polyfit*/
    protected int K = numBatch / 4;

    // why numbatch/4 ? --> the signal spectral is computed using as the mean of two samples
    // therefore:
    // divide for 2, because having numbatch batches, we consider numbatch / 2
    // divide again for 2, because we want to use only the first half of the spectr
    // (the second half is sigmetric)

    /**the order of the polynomial*/
    protected int polyOrder = 2;
    
    /**tells if estimated confidence interval is okay*/
    protected boolean confIntervalOk = false;
    protected boolean disableStatisticStop = false;
// ---------------------------------------------------------------------------------------------

// ---- NULL TEST HYPOTHESIS--------------------------------------------------------------------
    //@author Stefano Omini

    //NULL TEST is used to identify and stop all measures that are
    //very close to 0.0 (such measures cause troubles to the analyzer)

    //true if the analyzer must control whether the mean value is 0.0
    private boolean nullTestEnabled = true;
    //true if null test has determined that this measure is equal to zero
    private boolean measureIsZero = false;

    //contains the mean and the variance computed during the test
    private WeightedMeanVar nullTestData = new WeightedMeanVar();

    //the quantile required for the confidence interval of nullTest
    private double nullTestAlfa = 0.005;

    //to save time and resources, null test is not made for each sample, but repeated
    //every nullTestPeriod samples (testRate is equal to nullTestPeriod/maxData)
    private double nullTestRate = 0.01;
    //to save time and resources, null test is repeated every nullTestPeriod samples
    private int nullTestPeriod;

    //If measure is zero, this is the value of the upper limit of confidence interval
    private static double nullMeasure_upperLimit = 0.0;
// ---------------------------------------------------------------------------------------------

    /**
     *  Creates a new Dynamic Data Analyzer
     *  @param alpha the quantile required for the confidence interval
     *  @param precision maximum relative error
     *  @param maxData maximum number of data to be analyzed (at least 5000)
     */
    public NewDynamicDataAnalyzer(double alpha, double precision, int maxData) {
        if (alpha <= 0 || alpha >= 1 || precision <= 0 || precision >= 1 || maxData < 5000)
            throw new IllegalArgumentException("Wrong  Dynamic analyzer parameter  " + alpha + "  " + precision + "  " + maxData);
        this.alpha = alpha;
        this.precision = precision;
        this.maxData = maxData;

        // Initialize variables
        hMeans = new double[checkSize * hArrayMult];

        //@author Stefano Omini
        if (nullTestEnabled) {
            //initializes nullTestPeriod (the period of null test repetition)
            nullTestPeriod = (int) (nullTestRate * maxData);
        }
    }



    /**
     * Adds the new sample to the statistic.
     *
     * @param newSample the new sample
     * @param weight    the weight of the newSample, if it is not needed put 1.
     * @return true if the computation of confidence interval has finished (i.e. if the
     *         confidence interval is smaller than the one required by
     *         the user, or if the analyzed data are too many), false otherwise
     */
    public boolean addSample(double newSample, double weight) {
        // Skips 0 weight samples
        if (weight == 0.0)
            return end;

        // Abort simulation when maxData is reached
        if (nsamples > maxData && !confIntervalOk) {
            success = false;
            end = true;
            // Free allocated memory
            hMeans = null;
            return true;
        }
        // Max samples are reached and confidence interval too, but statistic stop was disabled
        else if (nsamples > maxData && confIntervalOk) {
            success = true;
            end = true;
            return true;
        }

        // Controls if analysis wasn't already terminated
        if (end)
            return true;

        nsamples++;

        //@author Stefano Omini
        //Measures which are often equal to zero are problematic for the
        //dynamic analyzer: in fact, in such cases, Schruben test for stationarity
        //may have no success.
        //For this reason a null test can be done in order to recognize if
        //a measure is almost always 0.

        if (nullTestEnabled) {

            //updates mean and var
            nullTestData.putNewSample(newSample, weight);

            //nullTest must be repeated every nullTestPeriod samples
            if ((nsamples % nullTestPeriod == 0)) {
                if (nullTest()) {
                    //null test is true (the measure is almost always equal to 0.0)
                    end = true;
                    success = true;
                    measureIsZero = true;
                    return end;
                }
            }
        }

        if (!delet)
            //transient not deleted
            batches.add(newSample, weight);

        // Execute heuristic
        if (!heuristicPassed) {
            heuristicPassed = FishmanRule();
            if (!heuristicPassed)
                return false;
            else {
                // Free allocated memory
                hMeans = null;

                // Heuristic to find spectral analysis parameters
                minSamples = transientLen;

                numBatch = 1 << ((int) Math.ceil((Math.log(Math.sqrt(transientLen / 4))) / Math.log(2)));
			    if (numBatch < 64)
				    numBatch = 64;
			    batchLen = (int) Math.ceil(transientLen / ((double) numBatch * 4));

                // Next iteration will perform spectral analysis
                return false;
            }
        }
        // At this point Fisherman Rule was passed and biased data truncated. New samples will
        // be analyzed with spectral analisys method to find confidence interval

        if (!delet) {
            //transient detected but not deleted yet

            if (nsamples >= minSamples) {
                //min number of samples respected

                //Apply the heuristic rule MSR5 to discard the most biased samples of the
                //time series. It minimizes the Mean Square Error.
                if (deleteTransient()) {

                    //transient deleted
                    if ((batch == numBatch - 1) && (nsamples % batchLen == 0)) {

                        //Applies the spectral test to generate the Confidence Intervals.
                        //true if the precision has been reached
                        end = HWtest();

                    } else
                        end = false;

                } else {
                    //transient not deleted yet
                    end = false;
                }
            }

        } else {
            // Correct problems with last batch after transient removal
            if (batch >= batchMean.length) {
                batch = batchMean.length - 1;
                batchMean[batch] *= weightBatchMean[batch];
            }
            
            //transient detected and deleted
            batchMean[batch] += newSample * weight;
            weightBatchMean[batch] += weight;

            if (nsamples % batchLen == 0) {
                batchMean[batch] /= weightBatchMean[batch];

                if (batch == (numBatch - 1)) {
                    end = HWtest();
                } else
                    batch++;
            } else
                end = false;
        }
        return end;
    }

    /**
     * It uses the Fishman Means Rule
     * Inserts a new value of the sequence, and approximately detects the initial
     * transient period.
     * <br>See
     * <br>
     *  Fishman G.S. 1973 揅oncepts and Methods in discrete Event Digital Simulation