#include <iostream>
#include <iomanip>
#include <cstring>
#include <fstream>
#include <cmath>
#include <unistd.h>
#include <TFile.h>
#include <TTree.h>
#include <TRandom.h>
#include <TMath.h>
using namespace std;

void plot(){
     float x,y,z,vx,vy,vz,MBvar;
     string root_name;
     string enter;
     cout << "Give me the name of your .root file: ";
//     cin >> root_name;
root_name="0";

     if(root_name=='0')
     root_name="Maxwell.root";

     TFile *f = new TFile(root_name.c_str());

     	if(f->IsZombie())
        cerr << "File not exist!\n", exit(0);

//------------------Maxwell-Boltmann-----------------//

     TH1F *h4 = new TH1F("h4","Maxwell-Boltzmann",100,0,500000);
     TF1 *MB=new TF1("MB","(4/sqrt([0]))*pow(1.5,[1]/(2*[2]*[3]))*x*x*exp(-1*([1]/(2*[2]*[3]))*x*x)",0,500000);
     MB->SetParNames("Pi","Mass","Kb","Temp");
     MB->SetParameters(TMath::Pi(),3.35e-24,TMath::K()*1e+7,3e+02);
      for(int i=0; i<1000000; i++)
       {
        MBvar=MB->GetRandom();

        h4->Fill(MBvar);

      }
     TCanvas *c3=new TCanvas("c3","c3", 700, 900);
     c3->cd();
     h4->Draw();
     c3->Update();
     TPaveStats *stath4=(TPaveStats*) c3->GetPrimitive("stats");
     c3->Close();

//---------------------Data----------------------//

     f->ls();
     TTree *tree3= (TTree*)f->Get("tree"); 
     TH3F *h3 = new TH3F("h3","Radial distribution",10,-1,1,10,-1,1,10,-1,1);
     TH2F *h1 = new TH2F("h1","Position distribution",10,-1,1,10,-1,1);
     TH2F *h2 = new TH2F("h2","Velocity distribution",10,-1,1,10,-1,1);
     TCanvas* c0= new TCanvas("c0","c0", 700, 900);
     TCanvas* c1= new TCanvas("c1","c1", 700, 900);
     TCanvas* c2= new TCanvas("c2","c2", 700, 900);
     leg = new TLegend(0.0997067,0.9,0.283724,0.986765);
     c0->SetFillColor(17);
     c1->SetFillColor(17);
     c2->SetFillColor(17);


     c0->cd();
     c0->SetGrid();
     tree3->Draw("sqrt(x*x+y*y+z*z)>>h1"," sqrt(x*x+y*y+z*z)>0");
     h1->SetTitle("Radial distribution");
     h1->SetLineColor(2);
     h1->SetMarkerStyle(8);
     h1->SetMarkerColor(28);
     h1->SetMarkerSize(0.8);
     h1->SetFillColorAlpha(5, 0.70);
     h1->SetFillStyle(3010);
     h1->Draw("SAMEE");
     c0->Update();


     c1->cd();
     c1->SetGrid();
     tree3->Draw("sqrt(vx*vx+vy*vy+vz*vz)>>h2","sqrt(vx*vx+vy*vy+vz*vz)>0");
     h2->SetTitle("Speed distribution");
     h2->SetLineColor(kRed);
     h2->SetMarkerStyle(21);
     h2->SetMarkerColor(28);
     h2->SetMarkerSize(0.8);
     h2->SetFillColorAlpha(5, 0.70);
     h2->SetFillStyle(3021);
     h2->Draw("SAMEE");
     h4->Draw("SAME");
     leg->SetFillColor(20);
     leg->AddEntry(h2,"Data (h2)","lp");
     leg->AddEntry(h4,"Maxwell-Boltzmann (h4)","l");
     leg->Draw();
     stath4->SetY1NDC(0.614706);
     stath4->SetY2NDC(0.775);
     stath4->SetX1NDC(0.779326);
     stath4->SetX2NDC(0.980205);
     stath4->Draw();
     c1->Update();

     c2->cd();
     c2->SetGrid();
     gStyle->SetPalette(55);
     tree3->Draw("x:y:z:8*sqrt(x*x+y*y+z*z)>>h3"," sqrt(x*x+y*y+z*z)>0","colz",100000);
     h3->SetTitle("Representation of gas in the sphere");
     c2->Update();




  /*   cout << "Press any key to exit..." << endl;

        system("read");
        f->Close();
	c0->Close();
	c1->Close();
	c2->Close();
        exit(0);*/
       

}