void histogram () {

   TCanvas *c = new TCanvas("Toy study","Toy Study", 1000, 700);
   c->Divide(2,2);
	
   c->cd(1);

   TF1 *f1= new TF1("f1","gaus(x)",-100,100);

   double mu_gen=0;
   double sigma=10;
   double cons=100;

   TH1F *h=new TH1F("h", "My histogram ; Gaussian mean 0 ; Number of entries",100,-60,60);
   TH1F *h_pull=new TH1F("h_pull","Pull distribution; (mu_gen-mu_fit)/mu_error; Number of entries ", 100, -4,4);
   TH1F *h_mu_difference= new TH1F("h_mu_difference","Mu_gen-Mu_fit; mu_gen-mu_fit; Number of entries",100,-0.1,0.1);
   TH1F *h_error= new TH1F("h_error","Error; error of mu; Number of entries ",100,-0.1,0.2);

   for (int j=0; j<10000; j++){
//std::cout<<"iterations "<<j<<std::endl;

      if (j>0) h->Reset();
      // set initial parameter in the function
      f1->SetParameter(0,cons);
      f1->SetParameter(1,mu_gen);
      f1->SetParameter(2,sigma);

      h->FillRandom("f1",10000);

      h->Fit("f1","R Q0 L");	
	
      double mu_fit=f1->GetParameter(1);
      double mu_error=f1->GetParError(1);
      double pull=( (mu_gen-mu_fit)/mu_error);

//std::cout<<"pull distribution"<<pull<<std::endl;

      double difference= mu_gen-mu_fit;
	
      h_mu_difference->Fill(difference);
      h_error->Fill(mu_error);
      h_pull->Fill(pull);

      f1->SetParameter(0,cons);
      f1->SetParameter(1,mu_gen);
      f1->SetParameter(2,sigma);



//TFitResultPtr r = h->Fit(f1);
//Int_t fitStatus = r;

   }
//std::cout<<"mu_fit"<<mu_fit<<std::endl;
   //double error=f1->GetParError(1);
   //legend-> AddEntry(	

	
   c->Draw();	
   h->Draw();

   c->cd(2);

   h_pull->Draw();

   c->cd(3);

   h_error->Draw();

   c->cd(4);

   h_mu_difference->Draw();

}