//macro to read mini root data
//and to book some histograms that will be stored in a root file
{

// Clear ROOT environment
        gROOT.Reset();
// Switch on ROOT timer
        gROOT.Time();
        gROOT->SetStyle("Plain");

        Int_t MaxEvents = 50000;
// Create a T49Run object
        T49Run *run = new T49Run();
        run->SetVerbose(kFALSE);


int how=0;

        Char_t *RunType="std+-160GeV-central-00B";
        Int_t StartAtRun = 1398;   //from 1389 - different trigger 
        Int_t EndAtRun = 1458;
        Int_t eventCount = 0;

// Define a pointer to an T49EventRoot object
        T49EventRoot    *event;
// Define a pointer to an T49ParticleRoot object
        T49ParticleRoot *track;


TH2F *hist1 = new TH2F("hist1", "-2.0<Ycms<-1.8", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist2 = new TH2F("hist2", "-1.8<Ycms<-1.6", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist3 = new TH2F("hist3", "-1.6<Ycms<-1.4", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist4 = new TH2F("hist4", "-1.4<Ycms<-1.2", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist5 = new TH2F("hist5", "-1.2<Ycms<-1.0", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist6 = new TH2F("hist6", "-1.0<Ycms<-0.8", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist7 = new TH2F("hist7", "-0.8<Ycms<-0.6", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist8 = new TH2F("hist8", "-0.6<Ycms<-0.4", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist9 = new TH2F("hist9", "-0.4<Ycms<-0.2", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist10 = new TH2F("hist10", "-0.2<Ycms<0.0", 100, -180, 180, 100, 0, 1.5);

TH2F *hist11 = new TH2F("hist11", "0.0<Ycms<0.2", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist12 = new TH2F("hist12", "0.2<Ycms<0.4", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist13 = new TH2F("hist13", "0.4<Ycms<0.6", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist14 = new TH2F("hist14", "0.6<Ycms<0.8", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist15 = new TH2F("hist15", "0.8<Ycms<1.0", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist16 = new TH2F("hist16", "1.0<Ycms<1.2", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist17 = new TH2F("hist17", "1.2<Ycms<1.4", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist18 = new TH2F("hist18", "1.4<Ycms<1.6", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist19 = new TH2F("hist19", "1.6<Ycms<1.8", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist20 = new TH2F("hist20", "1.8<Ycms<2.0", 100, -180, 180, 100, 0, 1.5 );

TH2F *hist21 = new TH2F("hist21", "2.0<Ycms<2.2", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist22 = new TH2F("hist22", "2.2<Ycms<2.4", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist23 = new TH2F("hist23", "2.4<Ycms<2.6", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist24 = new TH2F("hist24", "2.6<Ycms<2.8", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist25 = new TH2F("hist25", "2.8<Ycms<3.0", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist26 = new TH2F("hist26", "3.0<Ycms<3.2", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist27 = new TH2F("hist27", "3.2<Ycms<3.4", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist28 = new TH2F("hist28", "3.4<Ycms<3.6", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist29 = new TH2F("hist29", "3.6<Ycms<3.8", 100, -180, 180, 100, 0, 1.5 );
TH2F *hist30 = new TH2F("hist30", "3.8<Ycms<4.0", 100, -180, 180, 100, 0, 1.5 );



 //  open an output file
  TFile *hfile = new TFile("pt_phi_z_158_inY_negative.root", "RECREATE");    



// Loop over events in the data file
        while ((event = (T49EventRoot *) 
run->GetNextEvent(RunType,StartAtRun,EndAtRun)) && (eventCount < 
MaxEvents)) {

// find primary vertex
     T49Vertex *vertex = event->GetPrimaryVertex();

T49Vertex  *bpd = event->GetBPDVertex();

     double ntf=vertex->GetNTrkFit();
     double nto=vertex->GetNTrkOut();
     double ratio=ntf/nto;


//vertex cuts:
   if((vertex->GetPchi2()>0)
   &&(vertex->GetIflag()==0)
   &&(event->GetVertexX()>-0.1)
   &&(event->GetVertexX()<0.04)
   &&(event->GetVertexY()>-0.1)
   &&(event->GetVertexY()<0.13)
   &&(event->GetVertexZ()>-579.2)
   &&(event->GetVertexZ()<-578.6)
   &&(event->GetEveto()<10868.0)   //7.2 procent
   &&(nto!=0)
   &&(ratio>0.25))


{

how=how+1;
cout<<how<<endl;


//cout<<event->GetVertexZ()<<'\t'<<bpd->GetZ()<<endl;

// Get a list of all primary particles in the event
           TClonesArray *List = event->GetPrimaryParticles();
// Determine the number of primary particles
           Int_t max = List->GetEntries();


// Loop over all particles
           for(Int_t i=0;i<max;i++)
           {
               track = (T49ParticleRoot *) List->At(i);
 double np=track->GetNPoint();
 double nmp=track->GetNMaxPoint();

 double angle=(atan2(track->GetPy(), track->GetPx()))*(180.0/3.14);
 double rapidity=track->GetRap(0.1396);

 if(track->GetCharge()<0)
 {
   if(angle<0)
   angle=angle+360;
 angle=angle-180;  
 }


//track cuts:
    if(((track->GetIflag()&0xFF000000) == 0)
    &&(track->GetNMaxPoint()>30)
    &&(np/nmp>0.5)  
    &&(track->GetZFirst()<200.0)
    &&(fabs(track->GetBx())<2.0)
    &&(fabs(track->GetBy())<1.0)
    &&(track->GetCharge()<0)
    &&(track->GetPt()>0.005)
    &&(track->GetPt()<1.5))

   {

double ycms=rapidity-2.91;
double z=(rapidity-2.91)/2.91;



if(ycms>-2.0 && ycms<-1.8)
hist1->Fill(angle,track->GetPt());

if(ycms>-1.8 && ycms<-1.6)
hist2->Fill(angle,track->GetPt());

if(ycms>-1.6 && ycms<-1.4)
hist3->Fill(angle,track->GetPt());

if(ycms>-1.4 && ycms<-1.2)
hist4->Fill(angle,track->GetPt());

if(ycms>-1.2 && ycms<-1.0)
hist5->Fill(angle,track->GetPt());


if(ycms>-1.0 && ycms<-0.8)
hist6->Fill(angle,track->GetPt());

if(ycms>-0.8 && ycms<-0.6)
hist7->Fill(angle,track->GetPt());

if(ycms>-0.6 && ycms<-0.4)
hist8->Fill(angle,track->GetPt());

if(ycms>-0.4 && ycms<-0.2)
hist9->Fill(angle,track->GetPt());

if(ycms>-0.2 && ycms<0.0)
hist10->Fill(angle,track->GetPt());


if(ycms>0.0 && ycms<0.2)
hist11->Fill(angle,track->GetPt());

if(ycms>0.2 && ycms<0.4)
hist12->Fill(angle,track->GetPt());

if(ycms>0.4 && ycms<0.6)
hist13->Fill(angle,track->GetPt());

if(ycms>0.6 && ycms<0.8)
hist14->Fill(angle,track->GetPt());

if(ycms>0.8 && ycms<1.0)
hist15->Fill(angle,track->GetPt());


if(ycms>1.0 && ycms<1.2)
hist16->Fill(angle,track->GetPt());

if(ycms>1.2 && ycms<1.4)
hist17->Fill(angle,track->GetPt());

if(ycms>1.4 && ycms<1.6)
hist18->Fill(angle,track->GetPt());

if(ycms>1.6 && ycms<1.8)
hist19->Fill(angle,track->GetPt());

if(ycms>1.8 && ycms<2.0)
hist20->Fill(angle,track->GetPt());


if(ycms>2.0 && ycms<2.2)
hist21->Fill(angle,track->GetPt());

if(ycms>2.2 && ycms<2.4)
hist22->Fill(angle,track->GetPt());

if(ycms>2.4 && ycms<2.6)
hist23->Fill(angle,track->GetPt());

if(ycms>2.6 && ycms<2.8)
hist24->Fill(angle,track->GetPt());

if(ycms>2.8 && ycms<3.0)
hist25->Fill(angle,track->GetPt());


if(ycms>3.0 && ycms<3.2)
hist26->Fill(angle,track->GetPt());

if(ycms>3.2 && ycms<3.4)
hist27->Fill(angle,track->GetPt());

if(ycms>3.4 && ycms<3.6)
hist28->Fill(angle,track->GetPt());

if(ycms>3.6 && ycms<3.8)
hist29->Fill(angle,track->GetPt());

if(ycms>3.8 && ycms<4.0)
hist30->Fill(angle,track->GetPt());




    }
           }  
               eventCount++;
//  cout<<"end of event" << endl; 

        } 
 
} 


hfile->cd();
hist1->Write();
hist2->Write();
hist3->Write();
hist4->Write();
hist5->Write();
hist6->Write();
hist7->Write();
hist8->Write();
hist9->Write();
hist10->Write();
hist11->Write();
hist12->Write();
hist13->Write();
hist14->Write();
hist15->Write();
hist16->Write();
hist17->Write();
hist18->Write();
hist19->Write();
hist20->Write();
hist21->Write();
hist22->Write();
hist23->Write();
hist24->Write();
hist25->Write();
hist26->Write();
hist27->Write();
hist28->Write();
hist29->Write();
hist30->Write();


  hfile->Write();
  hfile->Close();

}