/* *
 * Title:        

 * Description:  

 * Copyright:    Copyright (c) 

 * Company:      

 * @author
 * @version 1.0
 */


package umontreal.iro.lecuyer.stochprocess;
import umontreal.iro.lecuyer.rng.*;
import umontreal.iro.lecuyer.probdist.*;
import umontreal.iro.lecuyer.randvar.*;



/**
 * This class represents a geometric variance gamma process S(t)
 * (see). This stochastic process is defined by the
 * equation
 * 

 * 
 * 

 * S(t) = S(0) exp(μt + X(t;σ, ν, θ) + ωt),
 * 

 * where X is a variance gamma process and
 * 

 * 
 * 

 * ω = (1/ν) ln(1 - θν - σ2ν/2).
 * 

 * 
 */
public class GeometricVarianceGammaProcess extends StochasticProcess  {
    protected VarianceGammaProcess vargamma;
    protected double        theta,
                            nu,
                            mu,
                            sigma,
                            omega,
                            muPlusOmega;
    protected double[]      mudt;



   /**
    * Constructs a new GeometricVarianceGammaProcess with parameters
    * 
    * θ = theta, 
    * σ = sigma, 
    * ν = nu,
    * 
    * μ = mu and initial value 
    * S(t0) = s0.
    * The stream is  used to generate the {@link VarianceGammaProcess} object used to implement
    * X in.
    * 
    */
   public GeometricVarianceGammaProcess (double s0, double theta,
                                         double sigma, double nu,
                                         double mu, RandomStream stream)  {
        vargamma = new VarianceGammaProcess (0.0, theta, sigma, nu, stream);
        setParams (s0, theta, sigma, nu, mu);
    }


   /**
    * Constructs a new GeometricVarianceGammaProcess.
    * The parameters 
    * θ, σ, ν are set to the parameters of the
    *  {@link VarianceGammaProcess} vargamma. The parameter μ
    * is set to mu and the initial values 
    * S(t0) = s0.
    * 
    */
   public GeometricVarianceGammaProcess (double s0, double mu,
                                         VarianceGammaProcess vargamma)
  {
        this.vargamma = vargamma;
        setParams (s0, vargamma.getTheta (), vargamma.getSigma (),
                   vargamma.getNu (), mu);
    }


   public double nextObservation() {
        double nextX  = vargamma.nextObservation();
        observationIndex = vargamma.getCurrentObservationIndex();
        // Could be different than simply 'observationIndex++' because of the
        // possibility of Gamma/Brownian bridge
        observationCounter++;

        double s = x0 * Math.exp (muPlusOmega * (t[observationIndex] - t[0])
                                  + nextX);
        path[observationIndex] = s;
        return s;
    }

   public double[] generatePath() {
        double s = x0;
        resetStartProcess();
        double[] vgpath = vargamma.generatePath();
        for (int i = 0; i < d; i++) {
            s *= Math.exp (mudt[i] + vgpath[i+1] - vgpath[i]);
            path[i+1] = s;
        }
        observationIndex = d;
        observationCounter++;
        return path;
    }
    // allows the user to create a path by specifiying the uniform random numbers to be used
   public double[] generatePath (double[] uniform01) {
        double s = x0;
        resetStartProcess();

        double[] vgpath = vargamma.generatePath(uniform01);
        for (int i = 0; i < d; i++) {
            s *= Math.exp (mudt[i] + vgpath[i+1] - vgpath[i]);
            path[i+1] = s;
        }
        observationIndex = d;
        observationCounter++;
        return path;
    }


    // method not verified by JS...  old stuff
   public double getCurrentUpperBound()  {
        // Find index for last observation generated (chronologically)
        int j = 0; // By default, t0 !
        int i = observationCounter - 1;
        double tForIthObserv;
        while (i > 0) {
            tForIthObserv = t[observationIndexFromCounter[i]];
            if (tForIthObserv <= t[observationCounter] && tForIthObserv > t[j])
                j = i;
            i--;
        }

        // Calculate bound following recipe
        double u = 0.0;
        GammaProcess gpos = ((VarianceGammaProcessDiff) vargamma).getGpos();
        double[] gposPath = gpos.getPath();
        double deltaGpos = gposPath[observationIndex] - gposPath[j];
        double s = path[observationIndex];
        if (muPlusOmega < 0)
             u = s * Math.exp (deltaGpos);
        else u = s * Math.exp (muPlusOmega * (t[observationIndex] - t[j])
                               + deltaGpos);
        return u;
    }


   /**
    * Resets the GeometricaVarianceGammaProcess,
    * but also applies the resetStartProcess method to the
    * {@link VarianceGammaProcess} object used to generate this process.
    * 
    */
   public void resetStartProcess()  {
        observationIndex   = 0;
        observationCounter = 0;
        vargamma.resetStartProcess();
    }


   /**
    * Sets the parameters
    * 
    * S(t0) = s0, 
    * θ = theta, 
    * σ = sigma,
    * 
    * ν = nu and 
    * μ = mu of the process.
    * Warning: This method will recompute some quantities stored internally,
    * which may be slow if called repeatedly.
    * 
    */
   public void setParams (double s0, double theta, double sigma, double nu,
                          double mu)  {
        this.x0    = s0;
        this.theta = theta;
        this.sigma = sigma;
        this.nu    = nu;
        this.mu    = mu;
        if (observationTimesSet) init(); // Otherwise no need to.
    }


   /**
    * Returns the value of the parameter θ.
    * 
    */
   public double getTheta()  { return theta; }


   /**
    * Returns the value of the parameter μ.
    * 
    */
   public double getMu()  { return mu; }


   /**
    * Returns the value of the parameter ν.
    * 
    */
   public double getNu()  { return nu; }


   /**
    * Returns the value of the parameter σ.
    * 
    */
   public double getSigma()  { return sigma; }


   /**
    * Returns the value of the quantity ω defined in.
    * 
    */
   public double getOmega()  { return omega; }


   /**
    * Returns a reference to the variance gamma process X defined
    * in the constructor.
    * 
    */
   public VarianceGammaProcess getVarianceGammaProcess()  {
        return vargamma;
}


    protected void init() {
        super.init();
        if (1 <= theta*nu + sigma*sigma*nu / 2.0)
           throw new IllegalArgumentException ("theta*nu + sigma*sigma*nu / 2 >= 1");
        omega = Math.log (1 - theta*nu - sigma*sigma*nu / 2.0) / nu;
        muPlusOmega = mu + omega;

        if (observationTimesSet) {
            // Telling the variance gamma proc. about the observ. times
            vargamma.setObservationTimes (t, d);

            // We need to know in which order the observations are generated
            this.observationIndexFromCounter
                = vargamma.getArrayMappingCounterToIndex();

            mudt = new double[d];
            for (int i = 0; i < d; i++)
                mudt[i] = muPlusOmega * (t[i+1] - t[i]);
        }
    }


   public void setStream (RandomStream stream)  {
        vargamma.setStream(stream);
   }


   public RandomStream getStream () {
      return vargamma.getStream();
   }

}