#include "peg_pic.h"


/**** Reading "general" parameters from the input file (pargv[1]) *******/

void rgen(char *pargv[])
{
	int    i, j, ticol1, ticol2, ticol3, tnl1, tnl2, tnl3;
	int    tHdel1, tHdel2, tHdel3, tHdex1, tHdex2, tHdex3;
	int    tHdi1, tHdi2, tHdi3, thh1, thh2, thh3;
	int    tdh12, tdh13, tdh21, tdh23, tdh31, tdh32;
	float  tpr1, tpr2, tpr3, tgt1, tgt2, tgt3, teex1, teex2, teex3;
	float  tdl1, tdl2, tdl3, tsel1, tsel2, tsel3, tspex1, tspex2, tspex3;
	float  teel1, teel2, teel3, tspel1, tspel2, tspel3, tnm1, tnm2, tnm3;
	float  tepel1, tepel2, tepel3, tacx1, tacx2, tacx3;
	float  tbcx1, tbcx2, tbcx3, tasc1, tasc2, tasc3;
	float  tbsc1, tbsc2, tbsc3, tacx12, tacx13, tacx21, tacx23, tacx31, tacx32;
	float  tbcx12, tbcx13, tbcx21, tbcx23, tbcx31, tbcx32;
	float  tasc12, tasc13, tasc21, tasc23, tasc31, tasc32;
	float  tbsc12, tbsc13, tbsc21, tbsc23, tbsc31, tbsc32;
	float  tepex1, tepex2, tepex3, tsi1, tsi2, tsi3, tsex1, tsex2, tsex3;
	float  tspi1, tspi2, tspi3, tei1, tei2, tei3, tepi1, tepi2, tepi3;
	char   aachar[80];

	InputDeck = fopen(pargv[1], "r");

	if (!InputDeck) {
		puts("\nCan't find input file ");
		puts(pargv[1]);
		exit(1);
	}

	/* Read lines until we get to numbers */

	while (fscanf(InputDeck,"%d %d %e %e %e %e %e %e %e",
			&nsp, &nc, &nc2p, &dt, &length, &area, &epsilon, &bmag, &psi) <9) {
		fscanf(InputDeck, "%s", aachar);
	}
	
	while (fscanf(InputDeck, "%e %e %d %d %d %e %d",
			&rhoback, &backj, &nfft, &nsmoothing, &MAXLEN, &dshort, &reflux) <7) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %d", &nfft1, &HISTMAX) <2) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%e %e %e %e %e %d %c", &dde, &extr, &extl, &extc, &extq, &dcramped, &src) <7) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%e %e %e %e %e", &dcbias, &ramp,
			&acbias, &w0, &theta0) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	/*
	 * ********* FIELD PARAMETERS *********
	 */
	while (fscanf(InputDeck, "%d %e %e", &field, &E_y, &x_E) <3) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %d %d %d %d %d %d", &slow, &dwindow,
			&datstep, &dmpstep, &last_step, &origdmp, &origres) <7) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %d %d %d %e %e",  &mflag, &mvstep, &Ef_flag, &T_diag,
			&poth[0], &poth[1]) <6) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %d", &stopfl, &dn, &NDP) <3) {
		fscanf(InputDeck, "%s", aachar);
	}

	if (NDP) {
		xs_mid   = (float *)malloc(NDP*sizeof(float));
		xf_mid   = (float *)malloc(NDP*sizeof(float));

		for (j=0; j<NDP; j++) {
			xs_mid[j] = xf_mid[j] = 0.0;
		}
	}

	for (j=0; j<NDP; j++) {
		while (fscanf(InputDeck, "%e %e", &xs_mid[j], &xf_mid[j]) <2) {
			fscanf(InputDeck, "%s", aachar);
		}
		// printf("start= %f\n", xs_mid[j]);
	}


	while (fscanf(InputDeck, "%d %e %e %e %e %d %e %e", &sourcD, &Src, &L_s, &nc05, &ht, &tas, &tasr, &alpha_Six) <8) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %d %e %e %e %e %d %e %e", &sourcD1, &timed, &Src1, &L_s1, &nc051, &ht1, &tas1, &tasr1, &alpha1_Six) <8) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %d", &ccc, &n_av, &cc_coef, &cc_step) <4) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %d %d %d %d %d", &Nnum, &ecollisional,	&ticol1,  &ticol2,  &ticol3, &dshort1) <6) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%e %e %e %d %e", &tpr1,&tnm1,&tgt1,&tnl1,&tdl1) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%e %e %e %d %e", &tpr2,&tnm2,&tgt2,&tnl2,&tdl2) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%e %e %e %d %e", &tpr3,&tnm3,&tgt3,&tnl3,&tdl3) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdel1, &tsel1, &teel1,	&tspel1, &tepel1) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdex1, &tsex1, &teex1, &tspex1, &tepex1) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %d %e %e %e %e", &tHdi1, &efi, &tsi1, &tei1, &tspi1, &tepi1) <6) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdel2, &tsel2, &teel2, &tspel2, &tepel2) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdex2, &tsex2, &teex2, &tspex2, &tepex2) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdi2, &tsi2, &tei2, &tspi2, &tepi2) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdel3, &tsel3, &teel3, &tspel3, &tepel3) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdex3, &tsex3, &teex3, &tspex3, &tepex3) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tHdi3, &tsi3, &tei3, &tspi3, &tepi3) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &thh1, &tacx1, &tbcx1, &tasc1, &tbsc1) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tdh12, &tacx12, &tbcx12, &tasc12, &tbsc12) <5) {
		fscanf(InputDeck, "%s", aachar); 
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tdh13, &tacx13, &tbcx13, &tasc13, &tbsc13) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &thh2, &tacx2, &tbcx2, &tasc2, &tbsc2) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tdh21, &tacx21, &tbcx21, &tasc21, &tbsc21) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tdh23, &tacx23, &tbcx23, &tasc23, &tbsc23) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &thh3, &tacx3, &tbcx3, &tasc3, &tbsc3) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tdh31, &tacx31, &tbcx31, &tasc31, &tbsc31) <5) {
		fscanf(InputDeck, "%s", aachar);
	}

	while (fscanf(InputDeck, "%d %e %e %e %e", &tdh32, &tacx32, &tbcx32, &tasc32, &tbsc32) <5) {
		fscanf(InputDeck, "%s", aachar);
	}


	if (sourcD1 && (!sourcD)) {
		// For the second source we need the first one
		sourcD= 1;
		Src= 0.;
	}

	if (Nnum) {
		/* Allocate memory for carged-neutral collision parameters */
		icolsp= (int *)malloc(Nnum*sizeof(int));
		icollisional= (int **)malloc(Nnum*sizeof(int *)); 
		n_neut_L= (int *)malloc(Nnum*sizeof(int));
		Hyd= (int *)malloc(Nnum*sizeof(int));
		Hydex= (int *)malloc(Nnum*sizeof(int));
		Hydion= (int *)malloc(Nnum*sizeof(int));
		HH= (int **)malloc(Nnum*sizeof(int *));
		Nneut0= (float *)malloc(Nnum*sizeof(float));
		Tneut0= (float *)malloc(Nnum*sizeof(float));
		dlabda= (float *)malloc(Nnum*sizeof(float));
		selsmax= (float *)malloc(Nnum*sizeof(float));
		elsengy0= (float *)malloc(Nnum*sizeof(float));
		elsengy1= (float *)malloc(Nnum*sizeof(float));
		elsengy2= (float *)malloc(Nnum*sizeof(float));
		sextmax= (float *)malloc(Nnum*sizeof(float));
		ext0= (float *)malloc(Nnum*sizeof(float));
		extengy1= (float *)malloc(Nnum*sizeof(float));
		extengy2= (float *)malloc(Nnum*sizeof(float));
		sionmax= (float *)malloc(Nnum*sizeof(float));
		ionengy0= (float *)malloc(Nnum*sizeof(float));
		ionengy1= (float *)malloc(Nnum*sizeof(float));
		ionengy2= (float *)malloc(Nnum*sizeof(float));
		mn= (float *)malloc(Nnum*sizeof(float));
		achrgx= (float **)malloc(Nnum*sizeof(float *));
		bchrgx= (float **)malloc(Nnum*sizeof(float *));
		ascat= (float **)malloc(Nnum*sizeof(float *));
		bscat= (float **)malloc(Nnum*sizeof(float *));

		for (i=0; i<Nnum; i++) {
			HH[i] = (int *)malloc(Nnum*sizeof(int));
			icollisional[i] = (int *)malloc(Nnum*sizeof(int));
			achrgx[i] = (float *)malloc(Nnum*sizeof(float));
			bchrgx[i] = (float *)malloc(Nnum*sizeof(float));
			ascat[i] = (float *)malloc(Nnum*sizeof(float));
			bscat[i] = (float *)malloc(Nnum*sizeof(float));
		}

		/* save parameters */
		icolsp[0] = ticol1;
		Nneut0[0] = tpr1;
		mn[0] = tnm1;
		Tneut0[0] = tgt1;
		n_neut_L[0] = tnl1;
		dlabda[0] = tdl1;
		Hyd[0] = tHdel1;
		selsmax[0] = tsel1;
		elsengy0[0] = teel1;
		elsengy1[0] = tspel1;
		elsengy2[0] = tepel1;
		Hydex[0] = tHdex1;
		sextmax[0] = tsex1;
		ext0[0] = teex1;
		extengy1[0] = tspex1;
		extengy2[0] = tepex1;
		Hydion[0] = tHdi1;
		sionmax[0] = tsi1;
		ionengy0[0] = tei1;
		ionengy1[0] = tspi1;
		ionengy2[0] = tepi1;
		HH[0][0] = thh1;
		achrgx[0][0] = tacx1;
		bchrgx[0][0] = tbcx1;
		ascat[0][0] = tasc1;
		bscat[0][0] = tbsc1;

		if (Nnum > 1) {
			icolsp[1] = ticol2;
			Nneut0[1] = tpr2;
			mn[1] = tnm2;
			Tneut0[1] = tgt2;
			n_neut_L[1] = tnl2;
			dlabda[1] = tdl2;
			Hyd[1] = tHdel2;
			selsmax[1] = tsel2;
			elsengy0[1] = teel2;
			elsengy1[1] = tspel2;
			elsengy2[1] = tepel2;
			Hydex[1] = tHdex2;
			sextmax[1] = tsex2;
			ext0[1] = teex2;
			extengy1[1] = tspex2;
			extengy2[1] = tepex2;
			Hydion[1] = tHdi2;
			sionmax[1] = tsi2;
			ionengy0[1] = tei2;
			ionengy1[1] = tspi2;
			ionengy2[1] = tepi2;
			HH[1][0] = tdh12;
			achrgx[1][0] = tacx12;
			bchrgx[1][0] = tbcx12;
			ascat[1][0] = tasc12;
			bscat[1][0] = tbsc12;
			HH[0][1] = thh2;
			HH[1][1] = tdh21;
			achrgx[0][1] = tacx2;
			achrgx[1][1] = tacx21;
			bchrgx[0][1] = tbcx2;
			bchrgx[1][1] = tbcx21;
			ascat[0][1] = tasc2;
			ascat[1][1] = tasc21;
			bscat[0][1] = tbsc2;
			bscat[1][1] = tbsc21;

			if (Nnum == 3) {
				icolsp[2] = ticol3;
				Nneut0[2] = tpr3;
				mn[2] = tnm3;
				Tneut0[2] = tgt3;
				n_neut_L[2] = tnl3;
				dlabda[2] = tdl3;
				Hyd[2] = tHdel3;
				selsmax[2] = tsel3;
				elsengy0[2] = teel3;
				elsengy1[2] = tspel3;
				elsengy2[2] = tepel3;
				Hydex[2] = tHdex3;
				sextmax[2] = tsex3;
				ext0[2] = teex3;
				extengy1[2] = tspex3;
				extengy2[2] = tepex3;
				Hydion[2] = tHdi3;
				sionmax[2] = tsi3;
				ionengy0[2] = tei3;
				ionengy1[2] = tspi3;
				ionengy2[2] = tepi3;
				HH[2][0] = tdh13;
				achrgx[2][0] = tacx13;
				bchrgx[2][0] = tbcx13;
				ascat[2][0] = tasc13;
				bscat[2][0] = tbsc13;
				HH[2][1] = tdh23;
				achrgx[2][1] = tacx23;
				bchrgx[2][1] = tbcx23;
				ascat[2][1] = tasc23;
				bscat[2][1] = tbsc23;
				HH[0][2] = thh3;
				HH[1][2] = tdh31;
				HH[2][2] = tdh32;
				achrgx[0][2] = tacx3;
				achrgx[1][2] = tacx31;
				achrgx[2][2] = tacx32;
				bchrgx[0][2] = tbcx3;
				bchrgx[1][2] = tbcx31;
				bchrgx[2][2] = tbcx32;
				ascat[0][2] = tasc3;
				ascat[1][2] = tasc31;
				ascat[2][2] = tasc32;
				bscat[0][2] = tbsc3;
				bscat[1][2] = tbsc31;
				bscat[2][2] = tbsc32;
			}
		}
	}

} /* End of RGEN */

/***************************************************************/
