// doscatter.r -- Run scatter_demo1 a bunch of times on the same potential
//                with different energy values.

// argument list -- hwidth -- half-width of potential
//                  height -- max value of potential (pos or neg)
//                  type -- 1 ==> step; 2 ==> square; 3 ==> gaussian
//                  energymin,denergy,energymax -- specify energy range
//                         and interval

// graphics -- makes a plot of reflection and transmission probablilities
//             as a function of the specified energy range

// output -- returns a list with the transmission and reflection
//           probabilities plus ancillary information

// other routines -- This routine needs to have scatter1_demo.r loaded
//                   for it to work.

// suggested starting call -- b = doscatter(1,1,2,.5,.05,2);

doscatter = function(hwidth,height,type,energymin,denergy,energymax) {
	if (nargs != 6) {
		error("Usage: doscatter(hwidth,height,type,energymin,denergy,energymax)");
	}

// set the integration step and the number of steps
	nsteps = 1001;
	dx = 0.02;

// set the energy values used
	num = int((energymax - energymin)/denergy + 1.5);
	range = [1:num]';
	erange = (range - 1)*denergy + energymin;

// initialize the vectors containing the reflection and transmission fractions
	refl = range;
	trans = range;

// run the loop
	for (i in 1:num) {
		energy = energymin + denergy*i;
		a = scatter1(nsteps,dx,height,hwidth,type,energy,0);
		refl[i] = a.reflprob;
		trans[i] = a.tranprob;
		printf(".");
	}
	printf("\n");

// plot the information
	pgopen();
	pgclear();
	pgwindow(1,0,0,15,10);
	pgaxes(1,0,erange[num],11,"energy",0,1,11,"trans(red) refl(blue) fracts");
	pgline(1,5,erange,trans);
	pgline(1,13,erange,refl);
	sprintf(xxx,"Potential: hwidth = %6.2f; height = %6.2f; type = %d",hwidth,height,type);
	pglabel(1,.1,.8,xxx);
	pgpause();
	pgclose();

// return the information
	return <<erange = erange; refl = refl; trans = trans; height = height; hwidth = hwidth; type = type>>;
}