{
	TFile* f = new TFile("sim_digi.root");
   	TTree* t = (TTree*)f->Get("er");
   	//t->AddFriend("er","reco.root");

   	TCanvas *c1 = new TCanvas("c1", "Telescope simulation", 1400, 1000);
   	c1->Divide(2,2);
   	c1->cd(1);
   	t->Draw("ERQTelescopeSiPoint_T1_DoubleSi_SD2_XY_0_X.fPDG","ERQTelescopeSiPoint_T1_DoubleSi_SD2_XY_0_X.fPDG == 1000020030");
   	c1->cd(2);
   	t->Draw("ERQTelescopeSiPoint_T1_DoubleSi_SD2_XY_0_X.fPDG","ERQTelescopeSiPoint_T1_DoubleSi_SD2_XY_0_X.fPDG == 1000010030");
   	c1->cd(3);
   	t->Draw("ERQTelescopeSiPoint_T2_DoubleSi_SD2_XY_1_X.fPDG","ERQTelescopeSiPoint_T2_DoubleSi_SD2_XY_1_X.fPDG == 1000020030");
   	c1->cd(4);
   	t->Draw("ERQTelescopeSiPoint_T2_DoubleSi_SD2_XY_1_X.fPDG","ERQTelescopeSiPoint_T2_DoubleSi_SD2_XY_1_X.fPDG == 1000010030");

	TCanvas *c2 = new TCanvas("c2", "Telescope de/E on digi info", 1400, 1000);
	c2->Divide(2,1);
	
	c2->cd(1);
	TH2F *T1deE = new TH2F("T1deE","T1 de(SingleSi) vs E(SingleSi+DoubleSi+CsI)",1000,0,0.25,1000,0,0.0015);
	T1deE->GetXaxis()->SetTitle("E");
	T1deE->GetYaxis()->SetTitle("de");
	t->Draw("ERQTelescopeSiDigi_T1_SingleSi_1_X_0.fEdep:ERQTelescopeSiDigi_T1_SingleSi_1_X_0.fEdep+(ERQTelescopeSiDigi_T1_DoubleSi_SD2_XY_0_X.fEdep+ERQTelescopeSiDigi_T1_DoubleSi_SD2_XY_0_Y.fEdep)/2.+ERQTelescopeCsIDigi_T1_CsI_1_0.fEdep>>T1deE","","gof");
	T1deE->Draw();

	c2->cd(2);
	TH2F *T2deE = new TH2F("T2deE","T1 de(DoubleSi) vs E(DoubleSi+CsI)",1000,0,0.25,1000,0,0.04);
	t->Draw("(ERQTelescopeSiDigi_T2_DoubleSi_SD2_XY_1_X.fEdep+ERQTelescopeSiDigi_T2_DoubleSi_SD2_XY_1_Y.fEdep)/2.:(ERQTelescopeSiDigi_T2_DoubleSi_SD2_XY_1_X.fEdep+ERQTelescopeSiDigi_T2_DoubleSi_SD2_XY_1_Y.fEdep)/2.+ERQTelescopeCsIDigi_T2_CsI_1_1.fEdep>>T2deE","","gof");
	T2deE->GetXaxis()->SetTitle("E");
	T2deE->GetYaxis()->SetTitle("de");
	T2deE->Draw();
}