Contour TGraph

ROOT
1

Hello all,
I have the following (probably quite elemntary) question, for those who might be able to help:
I scan a theoretical model’s parameter space (suppose 2-dimensional) to calculate corresponding chisquares among the model’s prediction and experimental data. This is done through an external personnal Fortran code.
I store the results in a file under the form (par1, par2, chisquare) and then read the file with a small root macro to create ntuples.
What I want to do, is to plot iso-chisquare contours on a 2-dimensional plane. To say this better, draw a smooth line connecting all the points with a given chisquare (of course the scan is done over a grid, so I guess a surface should first be generated in the 3D(par1, par2, chisquare) space, then a cut at a certain level and do a 2D projection of the cut).
If possible, I would also like to color the inside of the contours.

I guess I should use TGraphs (at first 3D and then 2D for the projection), since I’m not dealing with variables that could correspond to some sort of histogram.

Could somenoe help? By the way, to my knowledge there are really many theorists who have the same question!

Thanks a lot,
Andreas

2

Hi Andreas,

If your points are on a mesh Xi,Yi,Zi, i=0,i<N you can use a TGraph2D

TGraph2D g(filename); g.Draw(option); //see description of options in class doc
otherwise you can create a Tree with

TTree T("T","example"); T.ReadFile(filename,branchdescriptor); //see doc of TTree::readFile T.Draw("par2:par1","chisquare","cont0 colz");
Rene

3

Hello again,
Thank you very much for the reply.
I fear that this is, however, not what I’m looking for (and doesn’t really work, I guess I’m doing sth wrong). Perhaps if I described better what I want to do from scratch (since I’m not at all experienced with C++, perhaps I have it somewhere wrong):

I have, say, 4 input files in plain text format, with 4 columns each.

  • file1: m0, mhf, mchi, sigmav
  • file2: m0, mhf, mchi, mcharg
  • file3: m0, mhf, mchi, relic
  • file4: m0, mhf, mchi, chisquare

I want to plot, on the (m0, mhf) plane, the following contours:

  • Contours that join the points where mchi has one particular value, say mchi0
  • A contour that join the points where mcharg reaches a critical lower value, mcharg0
  • Contrours showing on which regions relic has a value between some limits, relic1 and relic2
  • One contour joining the points with a given chisquare, say chisquare0

The m0, mhf variables are indeed, as you said, on a grid. The rest are results from a computation wrt these parameters, so they have several non-canonical values.

I hope this time it is clearer…

As I said, I guess I should use a TGraph or TGraph2D and then some multigraph to add the contour graphs. But, I seem to miss the practical implementation…
Sorry I put so many things in one mail, but every time I seem to find some sort of “solution” (though never a real one) to one of them, this seems to conflict with the rest!

Thank you very very very much in advance!

Cheers,
Andreas

4

Hello,

I would like to bring up this topic once more, since I have a very similar problem.

I guess the solution would be to make a TGraph2D and plot contours from it. However, I miss a similar option as TH2D::SetContour(,), where I can specify which contours I want to have drawn.

What would be the way to do this?

Thanks for any help,
Jens

5

Can you post a small (running) example showing what you are doing ?

6

Well, not really, because I don’t know yet how to do it. But I can explain in simple words what I need:

I have a function z(x,y) given by a large number of discrete points (x,y,z). x and y are on a grid, while z can take arbitrary values. Now I would like to draw the contour in the (x,y) plane, on which z(x,y) is equal to a certain value z0.

I first thought I could create a TGraph2D from all the points and then use the “CONT” draw option. But how can I tell which contour is to be plotted then? I just need the one which fulfills z=z0.

Any ideas?

Thanks a lot,

Jens

7

If the (x,y) points are on a regular grid I would recommend you to use a TH2 . It will be more efficient. Then you can draw it with all the options described here: root.cern.ch/root/html/THistPainter.html

8

Hmm, I am not sure I understand this correctly. The problem is that the z values are arbitrary float numbers, since they are really values of the function z(x,y).

If z could just take integer values, I could set the bin contents of the histogram to the corresponding value and I would have what I need. But with z being float numbers?

Or did I miss something?

Thanks a lot!
Jens

9

If you use a TH2F or TH2D the bin contents (z values) can be any floating point number.

Rene

10

You can fill a bin with a non integer weight:
TH2::Fill(Double_t x, Double_t y, Double_t w)

11

Hello all,
I bring up this topic once again since (unfortunately!) I don’t seem to be solving the issue…

I do the following:

  • I create an NTuple with 3 branches, say mchi, mh, z. I want to plot iso-z contours on the (mchi, mh) plane: lines connecting all (mchi, mh) points with the same z.
  • I fill the histogram:
    for(Int_t i =0; iGetEntries();++i)
    {
    nt->GetEntry(i);
    hist->Fill(mchi, mh, z);
    }
    where hist is a TH2F.
  • I define contour levels:
    Double_t ContLevels[2];
    ContLevels[0] = 0.;
    ContLevels[1] = 1.;

histMIN->SetContour(2,ContLevels);

  • I plot:
    hist->Draw(“cont1z list”);
    c1->Update();

Unfortunately, the contours plotted have nothing to do with what I would like to have…
I don’t know exactly what is being plotted this way, it’s certainly not iso-z contours however, since for example z never reaches 1, nevertheless when I ask for a countour
at 1 a significant blob is plotted! I’m obviously doing sth very wrong, could someone help a bit?

As a second step, I’d like to put this contour in a TGraph along with various other points, lines etc (create a multigraph, which I have done several times without a problem), but I don’t seem to be able to get the contours in order to plot them afterwards…

Any ideas?

Thank you very much,
Andreas

12

Could you post the output of hist->Draw(“text”)?
or better provide the shortest possible running script showing your problem?

Rene

13

Hello,
Here’s the smallest possible script I can write, including the ntuple creation. As you’ll notice, I have imposed a countour at level 1, and this gives some contour, except that there is no “ratio” value going up to 1!
As an attachment, I also include the file I read in order to create my plots.
Thanks!
Andreas

{
gROOT->Reset();
ifstream DATA;
Int_t i = 0;
Float_t mchi, mh, mA0, mHCh, sigmav, resEIN, resNFW, resNFWc, ratio;
Int_t nlines = 0;

//…Read data file and create ntuple
DATA.open(“RESLinesbr-om0093-0128-3bd-0-mh500-dc120.txt”);
TFile *f = new TFile(“Einasto.root”,“RECREATE”);
TNtuple *ntEIN = new TNtuple(“ntEIN”,“data from ascii file”,“mchi:mA0:ratio”);

while (1) {
DATA>>mchi>>mh>>mA0>>mHCh>>sigmav>>resEIN>>resNFW>>resNFWc;
if (!DATA.good()) break;
ratio = sigmav/resEIN;
ntEIN->Fill(mchi, mA0, ratio);
}
DATA.close();
f->Write();

//…Create Canvas
TCanvas *c1 = new TCanvas();

//…Create Multigraph
//TMultiGraph *mg = new TMultiGraph();

TH2F *hist = new TH2F(“Einasto”, “Einasto”, 40, 40., 80., 60, 40., 100.);

f->Get(“ntEIN”);

ntEIN->SetBranchAddress(“mchi”,&mchi);
ntEIN->SetBranchAddress(“mA0”,&mA0);
ntEIN->SetBranchAddress(“ratio”,&ratio);

for(Int_t i =0; iGetEntries();++i)
{
ntEIN->GetEntry(i);
hist->Fill(mchi, mA0, ratio);
}

Double_t ContLevels[2];
ContLevels[0] = 0.;
ContLevels[1] = 1.;

hist->SetContour(2,ContLevels);

hist->Draw(“cont1z list”);
c1->Update();
}
RESLinesbr-om0093-0128-3bd-0-mh500-dc120.txt (102 KB)

14

I do not understand your problem. The contour z=1 works perfectly on your TH2. Simply prints the bin values as I suggested.

Rene

15

Hmmm, the point is that when I set the contour to be plotted at 1, I expect that nothing should be drawn, exactly because there is no (mchi, mA0) combination (point) that has a “ratio” value larger or equal to 1.

So, there is something wrong: I am obviously not telling ROOT to do what I want it to do.

16

This is not true. Simply draw your histogram with the option “text” or print the bin contents.

Rene

17

I see what you are saying.
I think the problem lies actually in my input file:
if right under
ntEIN->Fill(mchi, mA0, ratio);
I add a line
if (ratio > 1.) {cout<<ratio<<endl;}
I get no values printed on the terminal.

The problem is probably that there are more than 1 points in each bin (I am actually projecting a 4D parameter space on a 2D, so different points in 4D are the same in the 2D projection), so I guess that the total “ratio” value in each bin is the sum of the “ratio” contributions from each point in the bin, so in some cases this goes beyond 1.
I assume that the way to elliminate the contour would be to further divide the bin content by the number of points falling into the same bin.

Anyhow, now I would like to be able to retrieve this contour and draw it in a multigraph along with other things (for example, just scatter the points of my parameter space and then draw the contour on top of that, but also other things which actually require a TGraph of some sort), would that be possible?

Thank you very much for the help,
Andreas

18

Hi Andreas,

as far as I understand, you can use something like

TObjArray *contours = static_cast<TObjArray*>(gROOT->GetListOfSpecials()->FindObject("contours"));
Each entry in “contours” corresponds to one contour value and will be a TList which, in turn, contains at least one graph holding the points of a contour. There might be several graphs in the list if there is more than one contour at a value (i.e. several unconnected regions in the parameter space where your function takes the requested value).

Best,

Jens