/* This script Display the "dE vs E" pattern for HAND. The role of dE plane play scintillator in the preceding layer. The algoritm for choosing the correct bar is based on optical line from the target to the center of each bar.

*/




{


gROOT->Reset();






Double_t dist=6;	//distance in meters from the target

TFile *file0 = TFile::Open("~/E01015/Rootfiles/e01015_ND_1704.root");


//====================================== Geometrical structure of HAND =========================================================//


//Double_t xsizes[7]={0.11, 0.1, 0.125, 0.15, 0.25, 0.25, 0.25};// height (bar depth normal to beam direction) ---- ORIGINAL, without gaps!!!
Double_t xsizes[7]={0.112, 0.107, 0.13, 0.155, 0.252, 0.252, 0.252};// height (bar depth normal to beam direction)
Double_t zsizes[7]={0.02, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1};// Thickness (bar length parallel to beam direction)

int pads[7]={32, 30, 24, 22, 12, 12, 12};	//number of bars in each layer pads[0]=veto


Double_t xpos[7][32]={{-1.5646, -1.4526, -1.3406, -1.2286, -1.1166, -1.0046, -0.8926, -0.7806, -0.6686, -0.5566, -0.6115, -0.4995, -0.3875, -0.2755, -0.1635, -0.0515, 0.0606, 0.1726, 0.2846, 0.3966, 0.5086, 0.6206, 0.5646, 0.6766, 0.7886, 0.9006, 1.0126, 1.1246, 1.2366, 1.3486, 1.4606,1.5726},																{-1.5396, -1.4338, -1.328, -1.222, -1.1164, -1.0106, -0.9048, -0.799, -0.6932, -0.5874, -0.4816, -0.3758, -0.27, -0.1642,-0.0584, 0.0474, 0.1532, 0.259, 0.3648, 0.4706, 0.5764, 0.6822, 0.788, 0.8938, 0.9996, 1.1054, 1.2112, 1.317, 1.4228, 1.5286},{-1.5368, -1.406, -1.2752, -1.1444, -1.0136, -0.8828, -0.752, -0.6212, -0.4904, -0.3596, -0.2288, -0.0979, 0.0329, 0.1637, 0.2945, 0.4253,0.5561,0.6869,0.8177,0.9485,1.0793,1.2101,1.3409,1.4717},										{-1.4906, -1.3348, -1.179, -1.0232, -0.8674, -0.7116, -0.5683,-0.4375, -0.3067, -0.1759, -0.0576, 0.0482, 0.1665, 0.2973, 0.4281, 0.5589, 0.7022, 0.858, 1.0138, 1.1696, 1.3254, 1.4812},												{-1.4513, -1.192, -0.9328, -0.6736, -0.4143, -0.1551, 0.1042, 0.3634, 0.6226, 0.8819, 1.1411, 1.4004},					{-1.4513, -1.192, -0.9328, -0.6736, -0.4143, -0.1551, 0.1042, 0.3634, 0.6226, 0.8819, 1.1411, 1.4004},					{-1.4513, -1.192, -0.9328, -0.6736, -0.4143, -0.1551, 0.1042, 0.3634, 0.6226, 0.8819, 1.1411, 1.4004}};


Double_t zpos[7][32]={{-0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.25, -0.25, -0.25, -0.25, -0.25, -0.25, -0.25, -0.25, -0.25, -0.25, -0.25, -0.25, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15, -0.15,-0.15},                           {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},						{0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631, 0.1631},										{0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262, 0.3262},												{0.4893, 0.4893, 0.4893, 0.4893, 0.4893, 0.4893, 0.4893, 0.4893, 0.4893, 0.4893, 0.4893, 0.4893},					{0.6524, 0.6524, 0.6524, 0.6524, 0.6524, 0.6524, 0.6524, 0.6524, 0.6524, 0.6524, 0.6524, 0.6524},					{0.8155, 0.8155, 0.8155, 0.8155, 0.8155, 0.8155, 0.8155, 0.8155, 0.8155, 0.8155, 0.8155, 0.8155}};


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

TCanvas *c1 = new TCanvas("c1", "c1",15,49,1350,650);
c1.Divide(2,1);

int l;
Double_t h,x,angle;
char *line1,*line2,*YTitle,*XTitle,*cut,*dummy;
line1=(char*)malloc(1000);
line2=(char*)malloc(1000);
YTitle=(char*)malloc(1000);
XTitle=(char*)malloc(1000);
cut=(char*)malloc(1000);
dummy=(char*)malloc(1000);

sprintf(cut,"fEvtHdr.fEvtType==7");

TH2D *hr=new TH2D("hr","hr",100,0,4000,100,0,4000);
TH2D *hl=new TH2D("hl","hl",100,0,4000,100,0,4000);

T->SetMarkerStyle(3);


for(int i=1;i<5;i++)
	{
	cout<< "plane:" << i<<endl;
	for(int j=0;j<pads[i];j++)
		{
		angle=atan(xpos[i][j]/dist);	//particle angle
		x=TMath::Abs(tan(angle)*(zpos[i][j]-zpos[i-1][j]));		//displacement of xpos in the position of previous layer
			
		if(angle>=0) {h=xpos[i][j]-x;}
		else	{h=xpos[i][j]+x;}


		for(int k=0;k<pads[i-1];k++)


			{


			if(i==4)
				{

				if(k<=5 || k>=16 )
					{
					if(xpos[i-1][k]-xsizes[3]/2<h && h<xpos[i-1][k]+xsizes[3]/2)
						{
						sprintf(line1,"NA.nd.p%d.ra[%d]:NA.nd.p%d.ra[%d]>>hr",i-1,k,i,j);
						sprintf(line2,"NA.nd.p%d.la[%d]:NA.nd.p%d.la[%d]>>hl",i-1,k,i,j);

						cout << "line1:" << line1 << endl;

						sprintf(YTitle,"Layer %d, bar %d",i-1,k);
						sprintf(XTitle,"Layer %d, bar %d",i,j);

						c1.cd(1);
						T->Draw(line1,cut);
						

						hr->SetYTitle(YTitle);
						hr->SetXTitle(XTitle);
						c1.Paint();
						c1.Update();

						c1.cd(2);
						T->Draw(line2,cut);
						hl->SetYTitle(YTitle);
						hl->SetXTitle(XTitle);
						c1.Paint();
						c1.Update();
						gets(dummy);

						cout<<"Layer = "<<i<<" Bar = "<<j<<"||  Layer = "<<i-1<<"  Bar = "<<k<<endl;
						break;
						}
					}

				if((k>=6 && k<10) || (k>11 && k<16) )
					{
					if(xpos[i-1][k]-xsizes[2]/2<h && h<xpos[i-1][k]+xsizes[2]/2)
						{

						sprintf(line1,"NA.nd.p%d.ra[%d]:NA.nd.p%d.ra[%d]>>hr",i-1,k,i,j);
						sprintf(line2,"NA.nd.p%d.la[%d]:NA.nd.p%d.la[%d]>>hl",i-1,k,i,j);
						
						sprintf(YTitle,"Layer %d, bar %d",i-1,k);
						sprintf(XTitle,"Layer %d, bar %d",i,j);

						c1.cd(1);
						T->Draw(line1,cut);
						hr->SetYTitle(YTitle);
						hr->SetXTitle(XTitle);
						c1.Paint();
						c1.Update();

						c1.cd(2);
						T->Draw(line2,cut);
						hl->SetYTitle(YTitle);
						hl->SetXTitle(XTitle);
						c1.Paint();
						c1.Update();
						gets(dummy);

						cout<<"Layer = "<<i<<" Bar = "<<j<<"||  Layer = "<<i-1<<"  Bar = "<<k<<endl;
						break;
						}
					}

				if(k==10 || k==11)
					{
					if(xpos[i-1][k]-xsizes[1]/2<h && h<xpos[i-1][k]+xsizes[1]/2)
						{

						sprintf(line1,"NA.nd.p%d.ra[%d]:NA.nd.p%d.ra[%d]>>hr",i-1,k,i,j);
						sprintf(line2,"NA.nd.p%d.la[%d]:NA.nd.p%d.la[%d]>>hl",i-1,k,i,j);
						
						sprintf(YTitle,"Layer %d, bar %d",i-1,k);
						sprintf(XTitle,"Layer %d, bar %d",i,j);

						c1.cd(1);
						T->Draw(line1,cut);

						hr->SetYTitle(YTitle);
						hr->SetXTitle(XTitle);
						c1.Paint();
						c1.Update();

						c1.cd(2);
						T->Draw(line2,cut);
						hl->SetYTitle(YTitle);
						hl->SetXTitle(XTitle);
						c1.Paint();
						c1.Update();
						gets(dummy);

						cout<<"Layer = "<<i<<" Bar = "<<j<<"||  Layer = "<<i-1<<"  Bar = "<<k<<endl;
						break;
						}
					}
				
				} //if(i==4)


			else
				{
				if(xpos[i-1][k]-xsizes[i-1]/2<h && h<xpos[i-1][k]+xsizes[i-1]/2)
					{
					if(i==1)
						{
						sprintf(line1,"NA.veto.ra[%d]:NA.nd.p%d.ra[%d]>>hr",k,i,j);
						sprintf(line2,"NA.veto.la[%d]:NA.nd.p%d.la[%d]>>hl",k,i,j);
						}
					else
						{
						sprintf(line1,"NA.nd.p%d.ra[%d]:NA.nd.p%d.ra[%d]>>hr",i-1,k,i,j);
						sprintf(line2,"NA.nd.p%d.la[%d]:NA.nd.p%d.la[%d]>>hl",i-1,k,i,j);
						}
						
					sprintf(YTitle,"Layer %d, bar %d",i-1,k);
					sprintf(XTitle,"Layer %d, bar %d",i,j);
	
					c1.cd(1);
					T->Draw(line1,cut);
hl.SetMarkerStyle(2);
hr.SetMarkerStyle(2);
					hr->SetYTitle(YTitle);
					hr->SetXTitle(XTitle);
					c1.Paint();
					c1.Update();
	
					c1.cd(2);	
					T->Draw(line2,cut);

					hl->SetYTitle(YTitle);
					hl->SetXTitle(XTitle);
					c1.Paint();
					c1.Update();
					gets(dummy);

					cout<<"Layer = "<<i<<" Bar = "<<j<<"||  Layer = "<<i-1<<"  Bar = "<<k<<endl;
					break;
					}//if
				}
			}//for
		}//for
	}//for

}