ROOT logo

From $ROOTSYS/tutorials/hist/ContourList.C

// Getting Contours From TH2D
// Author: Josh de Bever
//         CSI Medical Physics Group
//         The University of Western Ontario
//         London, Ontario, Canada
//   Date: Oct. 22, 2004
//   Modified by O.Couet (Nov. 26, 2004)

TCanvas *ContourList(){

   const Double_t PI = TMath::Pi();

   TCanvas* c = new TCanvas("c","Contour List",0,0,600,600);
   c->SetRightMargin(0.15);
   c->SetTopMargin(0.15);

   Int_t i, j, TotalConts;

   Int_t nZsamples   = 80;
   Int_t nPhiSamples = 80;

   Double_t HofZwavelength = 4.0;       // 4 meters
   Double_t dZ             =  HofZwavelength/(Double_t)(nZsamples - 1);
   Double_t dPhi           = 2*PI/(Double_t)(nPhiSamples - 1);

   TArrayD z(nZsamples);
   TArrayD HofZ(nZsamples);
   TArrayD phi(nPhiSamples);
   TArrayD FofPhi(nPhiSamples);

   // Discretized Z and Phi Values
   for ( i = 0; i < nZsamples; i++) {
      z[i] = (i)*dZ - HofZwavelength/2.0;
      HofZ[i] = SawTooth(z[i], HofZwavelength);
   }

   for(Int_t i=0; i < nPhiSamples; i++){
      phi[i] = (i)*dPhi;
      FofPhi[i] = sin(phi[i]);
   }

   // Create Histogram
   TH2D *HistStreamFn = new TH2D("HstreamFn",
   "#splitline{Histogram with negative and positive contents. Six contours are defined.}{It is plotted with options CONT LIST to retrieve the contours points in TGraphs}",
   nZsamples, z[0], z[nZsamples-1], nPhiSamples, phi[0], phi[nPhiSamples-1]);

   // Load Histogram Data
   for (Int_t i = 0; i < nZsamples; i++) {
      for(Int_t j = 0; j < nPhiSamples; j++){
         HistStreamFn->SetBinContent(i,j, HofZ[i]*FofPhi[j]);
      }
   }

   gStyle->SetPalette(1);
   gStyle->SetOptStat(0);
   gStyle->SetTitleW(0.99);
   gStyle->SetTitleH(0.08);

   Double_t contours[6];
   contours[0] = -0.7;
   contours[1] = -0.5;
   contours[2] = -0.1;
   contours[3] =  0.1;
   contours[4] =  0.4;
   contours[5] =  0.8;

   HistStreamFn->SetContour(6, contours);

   // Draw contours as filled regions, and Save points
   HistStreamFn->Draw("CONT Z LIST");
   c->Update(); // Needed to force the plotting and retrieve the contours in TGraphs

   // Get Contours
   TObjArray *conts = (TObjArray*)gROOT->GetListOfSpecials()->FindObject("contours");
   TList* contLevel = NULL;
   TGraph* curv     = NULL;
   TGraph* gc       = NULL;

   Int_t nGraphs    = 0;
   Int_t TotalConts = 0;

   if (conts == NULL){
      printf("*** No Contours Were Extracted!\n");
      TotalConts = 0;
      return;
   } else {
      TotalConts = conts->GetSize();
   }

   printf("TotalConts = %d\n", TotalConts);

   for(i = 0; i < TotalConts; i++){
      contLevel = (TList*)conts->At(i);
      printf("Contour %d has %d Graphs\n", i, contLevel->GetSize());
      nGraphs += contLevel->GetSize();
   }

   nGraphs = 0;

   TCanvas* c1 = new TCanvas("c1","Contour List",610,0,600,600);
   c1->SetTopMargin(0.15);
   TH2F *hr = new TH2F("hr",
   "#splitline{Negative contours are returned first (highest to lowest). Positive contours are returned from}{lowest to highest. On this plot Negative contours are drawn in red and positive contours in blue.}",
   2, -2, 2, 2, 0, 6.5);

   hr->Draw();
   Double_t x0, y0, z0;
   TLatex l;
   l.SetTextSize(0.03);
   char val[20];

   for(i = 0; i < TotalConts; i++){
      contLevel = (TList*)conts->At(i);
      if (i<3) z0 = contours[2-i];
      else     z0 = contours[i];
      printf("Z-Level Passed in as:  Z = %f\n", z0);

      // Get first graph from list on curves on this level
      curv = (TGraph*)contLevel->First();
      for(j = 0; j < contLevel->GetSize(); j++){
         curv->GetPoint(0, x0, y0);
         if (z0<0) curv->SetLineColor(kRed);
         if (z0>0) curv->SetLineColor(kBlue);
         nGraphs ++;
         printf("\tGraph: %d  -- %d Elements\n", nGraphs,curv->GetN());

         // Draw clones of the graphs to avoid deletions in case the 1st
         // pad is redrawn.
         gc = (TGraph*)curv->Clone();
         gc->Draw("C");

         sprintf(val,"%g",z0);
         l.DrawLatex(x0,y0,val);
         curv = (TGraph*)contLevel->After(curv); // Get Next graph
      }
   }
   c1->Update();
   printf("\n\n\tExtracted %d Contours and %d Graphs \n", TotalConts, nGraphs );
   gStyle->SetTitleW(0.);
   gStyle->SetTitleH(0.);
   return c1;
}


Double_t SawTooth(Double_t x, Double_t WaveLen){

// This function is specific to a sawtooth function with period
// WaveLen, symmetric about x = 0, and with amplitude = 1. Each segment
// is 1/4 of the wavelength.
//
//           |
//      /\   |
//     /  \  |
//    /    \ |
//   /      \
//  /--------\--------/------------
//           |\      /
//           | \    /
//           |  \  /
//           |   \/
//

   Double_t y;
   if ( (x < -WaveLen/2) || (x > WaveLen/2)) y = -99999999; // Error X out of bounds
   if (x <= -WaveLen/4) {
      y = x + 2.0;
   } else if ((x > -WaveLen/4) && (x <= WaveLen/4)) {
      y = -x ;
   } else if (( x > WaveLen/4) && (x <= WaveLen/2)) {
      y = x - 2.0;
   }
   return y;
}
 ContourList.C:1
 ContourList.C:2
 ContourList.C:3
 ContourList.C:4
 ContourList.C:5
 ContourList.C:6
 ContourList.C:7
 ContourList.C:8
 ContourList.C:9
 ContourList.C:10
 ContourList.C:11
 ContourList.C:12
 ContourList.C:13
 ContourList.C:14
 ContourList.C:15
 ContourList.C:16
 ContourList.C:17
 ContourList.C:18
 ContourList.C:19
 ContourList.C:20
 ContourList.C:21
 ContourList.C:22
 ContourList.C:23
 ContourList.C:24
 ContourList.C:25
 ContourList.C:26
 ContourList.C:27
 ContourList.C:28
 ContourList.C:29
 ContourList.C:30
 ContourList.C:31
 ContourList.C:32
 ContourList.C:33
 ContourList.C:34
 ContourList.C:35
 ContourList.C:36
 ContourList.C:37
 ContourList.C:38
 ContourList.C:39
 ContourList.C:40
 ContourList.C:41
 ContourList.C:42
 ContourList.C:43
 ContourList.C:44
 ContourList.C:45
 ContourList.C:46
 ContourList.C:47
 ContourList.C:48
 ContourList.C:49
 ContourList.C:50
 ContourList.C:51
 ContourList.C:52
 ContourList.C:53
 ContourList.C:54
 ContourList.C:55
 ContourList.C:56
 ContourList.C:57
 ContourList.C:58
 ContourList.C:59
 ContourList.C:60
 ContourList.C:61
 ContourList.C:62
 ContourList.C:63
 ContourList.C:64
 ContourList.C:65
 ContourList.C:66
 ContourList.C:67
 ContourList.C:68
 ContourList.C:69
 ContourList.C:70
 ContourList.C:71
 ContourList.C:72
 ContourList.C:73
 ContourList.C:74
 ContourList.C:75
 ContourList.C:76
 ContourList.C:77
 ContourList.C:78
 ContourList.C:79
 ContourList.C:80
 ContourList.C:81
 ContourList.C:82
 ContourList.C:83
 ContourList.C:84
 ContourList.C:85
 ContourList.C:86
 ContourList.C:87
 ContourList.C:88
 ContourList.C:89
 ContourList.C:90
 ContourList.C:91
 ContourList.C:92
 ContourList.C:93
 ContourList.C:94
 ContourList.C:95
 ContourList.C:96
 ContourList.C:97
 ContourList.C:98
 ContourList.C:99
 ContourList.C:100
 ContourList.C:101
 ContourList.C:102
 ContourList.C:103
 ContourList.C:104
 ContourList.C:105
 ContourList.C:106
 ContourList.C:107
 ContourList.C:108
 ContourList.C:109
 ContourList.C:110
 ContourList.C:111
 ContourList.C:112
 ContourList.C:113
 ContourList.C:114
 ContourList.C:115
 ContourList.C:116
 ContourList.C:117
 ContourList.C:118
 ContourList.C:119
 ContourList.C:120
 ContourList.C:121
 ContourList.C:122
 ContourList.C:123
 ContourList.C:124
 ContourList.C:125
 ContourList.C:126
 ContourList.C:127
 ContourList.C:128
 ContourList.C:129
 ContourList.C:130
 ContourList.C:131
 ContourList.C:132
 ContourList.C:133
 ContourList.C:134
 ContourList.C:135
 ContourList.C:136
 ContourList.C:137
 ContourList.C:138
 ContourList.C:139
 ContourList.C:140
 ContourList.C:141
 ContourList.C:142
 ContourList.C:143
 ContourList.C:144
 ContourList.C:145
 ContourList.C:146
 ContourList.C:147
 ContourList.C:148
 ContourList.C:149
 ContourList.C:150
 ContourList.C:151
 ContourList.C:152
 ContourList.C:153
 ContourList.C:154
 ContourList.C:155
 ContourList.C:156
 ContourList.C:157
 ContourList.C:158
 ContourList.C:159
 ContourList.C:160
 ContourList.C:161
 ContourList.C:162
 ContourList.C:163
 ContourList.C:164
 ContourList.C:165
 ContourList.C:166
 ContourList.C:167
 ContourList.C:168
 ContourList.C:169
 ContourList.C:170
 ContourList.C:171
 ContourList.C:172
 ContourList.C:173
 ContourList.C:174