#include "Riostream.h"
#include "TH1.h"
#include "TCanvas.h"

void caen()
{
gROOT->Reset();
int Ext,i;
const int VectSize =256;
double W1[VectSize],W2[VectSize],W3[VectSize],I1[VectSize],W4[VectSize];
double ds1;
float f1,PSD,Qslow,Qtot,Qtot1,Qfast,P;
FILE *fpIn;

TCanvas *c1=new TCanvas("c1"," ",600,600);

TFile *f=new TFile("tree1.root","NEW");
TNtuple *ntuple=new TNtuple("ntuple"," ","Qtot:Qtot1:Qfast:Qslow");

TH1F *h1=new TH1F("h1","x distribution",VectSize,0,VectSize);
TH1F *h2=new TH1F("h2","x distribution",VectSize,0,VectSize);
TH1F *h3=new TH1F("h3","x distribution",VectSize,0,VectSize);
TH1F *h4=new TH1F("h4","x distribution",VectSize,0,VectSize);

c1->Divide(1,3);

Ext=0;
fpIn=fopen("24062012pu239sn113rec256thres2000mes17hour.txt","r");
if(fpIn==NULL) return;

do{ 
	for(i=0;i<VectSize;i++) {
			fscanf( fpIn, "%f\n",&f1 );
			W1[i]=f1;
			h1->SetBinContent(i,W1[i]);}
		    
		ds1=0;
		for(i=0;i<40;i++) {
		ds1+=W1[i];	
		};

		for(i=0;i<VectSize;i++) {
		W1[i]-=ds1/40.0;
		}								

		for(i=0;i<VectSize;i++) {W3[i]=(-1)*W1[i];
		h2->SetBinContent(i,W3[i]); }
				
		for(i=0;i<VectSize;i++) {if (W3[i]>0) W2[i]=W3[i];}  
		
		P=W3[0];
	        for(i=0;i<VectSize;i++) {
		if(P<W3[i]) P=W3[i];}
		
		for(i=0;i<VectSize;i++) { W4[i]=(W3[i]/P); h4->SetBinContent(i,W4[i]);}
		
		
			Qtot=0;
			for(i=40;i<250;i++) {
			Qtot+=W4[i];
			}
			
			Qfast=0;
			for(i=40;i<60;i++) {
			Qfast+=W4[i];}
			
			Qslow=0;
			for(i=60;i<70;i++) {
			Qslow+=W4[i];}

			PSD=(Qtot-Qfast)/Qtot;
			ntuple->Fill(Qtot,PSD,Qfast,Qslow);
		
			
printf("%d\n",Ext);
Ext++;		
}
while(!feof(fpIn));
fclose(fpIn);

c1->cd(1);
h1->Draw();
c1->cd(2);
h2->Draw();
c1->cd(3);
h4->Draw();
c1->Update();

f->Write();
ntuple->Write();
}