/*

This script calculate the desired voltages for HAND scintillators in two cases:

1. Cosmic run=1;
2. Beam (run=0)

The difference is the desired ADC channel at each layer. For cosmic, the energy deposit is different between layer because of different heights of bars. The ADC channel location is added based on the minimum ionizing particlas with 2.2 MeV/cm Energy deposit.

In the beam case, all bars have same width, however the kinetic energy of protons is decreasing as we moving from layer 1 to 6. Correct values for the ADC channels must be added to the code.


This script must be run with "read_gain.C"
and gain curve parameter must be available in the corresponding .dat file.


The output of this script are two files, one for left side and second for the right side, with voltage values.


*/






#include "read_gains.C"
#include <iostream>

void desired_voltage()	{

const int nbars[7]={32,30,24,22,12,12,12};



//======================== read_gains.C ====================
double  c_EL[24]={0},c_ER[24]={0},
	c_dEL[24]={0},c_dER[24]={0};	//  K*V**n+const  :   K = c_ , n = p_ , const = cnst_

double  cnst_EL[24]={0},cnst_ER[24]={0},
	cnst_dEL[24]={0},cnst_dER[24]={0};

int	p_EL[24],p_ER[24],
	p_dEL[24],p_dER[24];
//============================================================

double VER[24]={0}, VdER[24]={0};
double VEL[24]={0}, VdEL[24]={0};


double DEL,DER,DdEL,DdER, A;

double ADC_E=500;
double ADC_dE=250;


read_gains(	&c_EL[], &c_ER[], &c_dEL[], &c_dER[],
		&cnst_EL[], &cnst_ER[], &cnst_dEL[], &cnst_dER[],
		&p_EL[], &p_ER[], &p_dEL[], &p_dER[]);

//cout<<"cnst = "<< cnst_R4[6]<<endl;

for(int j=0;j<24;j++)
	{
	
	DEL=ADC_E-cnst_EL[j];
	DER=ADC_E-cnst_ER[j];
	DdEL=ADC_dE-cnst_dEL[j];
	DdER=ADC_dE-cnst_dER[j];

	if(DEL>=0 && c_EL[j]>0)	{VEL[j]=(DEL/c_EL[j])**(1./p_EL[j]);}
	else	cout<<" Problem with LEFT PMT "<<j<<endl;
	if(DER>=0 && c_ER[j]>0)	{VER[j]=(DER/c_ER[j])**(1/p_ER[j]);}
	else	cout<<" Problem with RIGHT PMT "<<j<<endl;	

	if(DdEL>=0 && c_dEL[j]>0)	{VdEL[j]=(DdEL/c_dEL[j])**(1./p_dEL[j]);}
	else	cout<<" Problem with LEFT PMT "<<j<<endl;
	if(DdER>=0 && c_dER[j]>0)	{VdER[j]=(DdER/c_dER[j])**(1/p_dER[j]);}
	else	cout<<" Problem with RIGHT PMT "<<j<<endl;	
	


	} // for j



//=================================================================================================

char lineL[200], lineR[200], valEL[200], valER[200],valdEL[200],valdER[200];

ofstream fileEL, fileER,filedEL,filedER;
remove("Voltages_Eleft.dat");
remove("Voltages_Eright.dat");
remove("Voltages_dEleft.dat");
remove("Voltages_dEright.dat");


fileEL.open("Voltages_Eleft.dat", ios::app);
fileER.open("Voltages_Eright.dat",ios::app);
filedEL.open("Voltages_dEleft.dat", ios::app);
filedER.open("Voltages_dEright.dat",ios::app);

for(int i=0;i<24;i++)
	{

	sprintf(valEL,"%.1f ",VEL[i]);
	sprintf(valER,"%.1f ",VER[i]);
	sprintf(valdEL,"%.1f ",VdEL[i]);
	sprintf(valdER,"%.1f ",VdER[i]);

	fileEL<<valEL;
	fileER<<valER;
	filedEL<<valdEL;
	filedER<<valdER;
	}
fileEL<<endl;
fileER<<endl;
filedEL.close();
filedER.close();
		
}//void