void write_histograms()
{
   // Create a Tree with a few branches of type histogram
   // 25000 entries are filled in the Tree
   // For each entry, the copy of 3 histograms is written
   // The data base will contain 75000 histograms.
   gBenchmark->Start("write_histograms");
   TFile f("ht.root", "recreate");
   auto T     = new TTree("T", "test");
   auto hpx   = new TH1F("hpx", "This is the px distribution", 100, -4, 4);
   auto hpxpy = new TH2F("hpxpy", "py vs px", 40, -4, 4, 40, -4, 4);
   auto hprof  = new TProfile("hprof", "Profile of pz versus px", 100, -4, 4, 0, 20);
   T->Branch("hpx", "TH1F", &hpx, 32000, 0);
   T->Branch("hpxpy", "TH2F", &hpxpy, 32000, 0);
   T->Branch("hprof", "TProfile", &hprof, 32000, 0);
   Float_t px, py, pz;
   for (Int_t i = 0; i < 25000; i++) {
      if (i % 1000 == 0)
         printf("at entry: %d\n", i);
      gRandom->Rannor(px, py);
      pz = px * px + py * py;
      hpx->Fill(px);
      hpxpy->Fill(px ,py);
      hprof->Fill(px, pz);
      T->Fill();
   }
   T->Print();
   f.Write();
   gBenchmark->Show("write_histograms");
}

void read_histogram1()
{
   // Connect Tree generated by htw and show histograms for entry 12345
   auto f = TFile::Open("ht.root");
   auto T = f->Get<TTree>("T");
   TH1F *hpx = nullptr;
   TH2F *hpxpy = nullptr;
   TProfile *hprof = nullptr;
   T->SetBranchAddress("hpx", &hpx);
   T->SetBranchAddress("hpxpy", &hpxpy);
   T->SetBranchAddress("hprof", &hprof);
   T->GetEntry(12345);
   auto c1 = new TCanvas("c1", "test", 10, 10, 600, 1000);
   c1->Divide(1, 3);
   c1->cd(1);
   hpx->Draw();
   c1->cd(2);
   hpxpy->Draw();
   c1->cd(3);
   hprof->Draw();
   c1->Print("htr1.png");
}

void read_histogram2()
{
   // Connect Tree generated by htw and show histograms for entry 12345
   // a variant of read_histogram1
   auto f = TFile::Open("ht.root");
   auto T = f->Get<TTree>("T");
   auto c1 = new TCanvas("c1", "test", 10, 10, 600, 1000);
   c1->Divide(1, 3);
   c1->cd(1);
   T->Draw("hpx.Draw()", "", "goff", 1, 12345);
   c1->cd(2);
   T->Draw("hpxpy.Draw()", "", "goff", 1, 12345);
   c1->cd(3);
   T->Draw("hprof.Draw()", "", "goff", 1, 12345);
   c1->Print("htr2.png");
}

void read_histogram3()
{
   // Connect Tree generated by htw
   // read all histograms and plot the RMS of hpx versus the Mean of hprof
   // for each of the 25000 entries
   auto f = TFile::Open("ht.root");
   auto T = f->Get<TTree>("T");
   auto c1 = new TCanvas("c1", "test", 10, 10, 600, 400);
   T->Draw("hpx.GetRMS():hprof.GetMean()");
   c1->Print("htr3.png");
}

void tree201_histograms()
{
   write_histograms();
   read_histogram1();
   read_histogram2();
   read_histogram3();
}