class WaveVector {
	int[] Re;
	int[] Im;
	int npts;

	public WaveVector() {
		this(256);
	}

	public WaveVector(int npts) {
		this.npts = npts;
		Re = new int[npts];
		Im = new int[npts];
		for(int i = 0; i<npts; i++) {
			Re[i] = 0;
			Im[i] = 0;
		}
	}


	double gaussian(double x, double mu, double sigma) {
		double d = (x - mu) / sigma;
		return Math.exp(-d*d/2) / sigma / Math.sqrt(2*Math.PI);
	}

	public void SetGaussian(double a0, double mu,
			 double xmin, double xmax, double v) {
		// Psi
		double l = Physics.h / Math.sqrt(2 * Physics.me * Physics.q * v);
		double a = l * a0;
		double k0 = 2*Math.PI / l;

		//estimate the max value of the wavevector, etc, based on initial amplitude.
		double scalef = ((double)(Integer.MAX_VALUE/4))
			/ Math.pow(2.0/(Math.PI*4*l*l), 0.25);

		//given initial params, setup the state
		double xrange = xmax - xmin;
		double dx = xrange/(npts-1);

		for(int i = 0; i<npts; i++) {
			double x = (i*dx + xmin) - mu;
			double mag = scalef
				* Math.exp(-(x*x)/(a*a))
				* Math.pow(2.0/(Math.PI * a * a), 0.25);
			Re[i] = (int)(mag * Math.cos(k0*x));
			Im[i] = (int)(mag * Math.sin(k0*x));
		}
	}

	public void init_wave(double sigma, double mu, double v, double dx) {
		double [] Psi_real = new double[npts];
		double [] Psi_imag = new double[npts];
		double [] probs = new double[npts];

		double lambda = Physics.h/(Physics.me*v);
		double k0 = 1 / (lambda * 2 * Math.PI);
		double tprob = 0;

		for(int i=0;i<npts;i++) {
			double x = (i - npts/2.0)*dx;
			double d = (x-mu)/sigma;
			double p = Math.exp(-0.5*d*d);
			probs[i] = p;
			tprob += p;
		}
		for(int i=0;i<npts;i++) {
			probs[i] /= tprob;
		}
		for(int i=0;i<npts;i++) {
			double x = (i - npts/2.0)*dx;
			double mag = Math.sqrt(probs[i]);
			Psi_real[i] = mag * Math.cos(k0*x);
			Psi_imag[i] = mag * Math.sin(k0*x);
		}

		double maxval = 0;
		for(int i=0;i<npts;i++) {
			maxval = Math.max(maxval, Math.abs(Psi_real[i]));
			maxval = Math.max(maxval, Math.abs(Psi_imag[i]));
		}
		/* We'll scale our integers so that they can hold values up to
		   almost twice the largest initial value before they incur
		   an overflow. */
		double scaleFactor = (1<<30)/maxval;

		/* Convert to integers. */
		for(int i=0;i<npts;i++) {
			Re[i] = (int)(Psi_real[i]*scaleFactor);
			Im[i] = (int)(Psi_imag[i]*scaleFactor);
		}
	}
}