using namespace std;

#include <cmath>
#include <cstdlib>
#include <iostream>
#include <fstream>

int main()
{
	// open files
	
	ifstream input("spin_08.paw");	// input file
	ofstream output("spin_dist08.paw");	// output file
	
	int l = 200;					// length of input & output file
	
	double energy[l], spin[l];	// vectors for input file. Units: energy: MeV
	double U[l], strength[l], strengtherr[l], sigma2[l];	// vectors for output file. Units: energy: MeV, RSF: MeV^(-3)
	
	const float a_lev = 10.33;	// const. factor in MeV^-1
	double A = 208;
	const float del = 2.226;
	
	// read from input to vectors
	
	int i = 0;
	while(input){
		input >> energy[i] >> spin[i];
		U[i] = energy[i]/1000.0 - del;
		sigma2[i] = 0.0888*pow(A, 0.667)*sqrt(a_lev*U[i]);
		strength[i] = ((2.0*spin[i]+1.0)/(2.0*sigma2[i]))*exp(-(spin[i]+0.5)*(spin[i] +0.5)/(2.00*sigma2[i]));
		output << energy[i]/1000.0 << '\t' << spin[i] << '\t' << strength[i]<< endl;
		i++;
	}
	
	input.close();
	output.close();
}