TGraph2D::GetHistogram leads to segfault

Hi

I have a rather simple code in pyROOT which leads to a segfault (but no backtrace) under certain circumstances. A minimal working example is found here (sorry for the length of the lists, but I could not reproduce the problem with shorter lists):

[code]#!/usr/bin/python

from ROOT import *
from array import array

gROOT.SetBatch(True)

gStyle.SetPaintTextFormat(‘5.2f’)

x=array(‘f’)
y=array(‘f’)
z=array(‘f’)

x.append(100.0)
x.append(100.0)
x.append(100.0)
x.append(112.5)
x.append(112.5)
x.append(117.5)
x.append(130.0)
x.append(132.5)
x.append(142.5)
x.append(150.0)
x.append(155.0)
x.append(157.5)
x.append(167.5)
x.append(175.0)
x.append(175.0)
x.append(182.5)
x.append(192.5)
x.append(200.0)
x.append(200.0)
x.append(200.0)
x.append(207.5)
x.append(217.5)
x.append(225.0)
x.append(225.0)
x.append(225.0)
x.append(232.5)
x.append(242.5)
x.append(250.0)
x.append(250.0)
x.append(250.0)
x.append(250.0)
x.append(257.5)
x.append(267.5)
x.append(275.0)
x.append(275.0)
x.append(275.0)
x.append(275.0)
x.append(282.5)
x.append(300.0)
x.append(300.0)
x.append(300.0)
x.append(300.0)
x.append(300.0)
x.append(325.0)
x.append(325.0)
x.append(325.0)
x.append(325.0)
x.append(330.0)
x.append(345.0)
x.append(350.0)
x.append(350.0)
x.append(350.0)
x.append(350.0)
x.append(362.5)
x.append(375.0)
x.append(375.0)
x.append(375.0)
x.append(387.5)
x.append(392.5)
x.append(400.0)
x.append(400.0)
x.append(400.0)
x.append(407.5)
x.append(412.5)
x.append(425.0)
x.append(425.0)
x.append(425.0)
x.append(437.5)
x.append(450.0)
x.append(450.0)
x.append(450.0)
x.append(462.5)
x.append(475.0)
x.append(475.0)
x.append(487.5)
x.append(500.0)
x.append(500.0)

y.append(0.0)
y.append(35.0)
y.append(65.0)
y.append(12.5)
y.append(47.5)
y.append(82.5)
y.append(30.0)
y.append(67.5)
y.append(107.5)
y.append(50.0)
y.append(50.0)
y.append(92.5)
y.append(132.5)
y.append(25.0)
y.append(75.0)
y.append(117.5)
y.append(157.5)
y.append(0.0)
y.append(100.0)
y.append(50.0)
y.append(142.5)
y.append(182.5)
y.append(125.0)
y.append(25.0)
y.append(75.0)
y.append(167.5)
y.append(207.5)
y.append(0.0)
y.append(100.0)
y.append(150.0)
y.append(50.0)
y.append(192.5)
y.append(232.5)
y.append(125.0)
y.append(175.0)
y.append(25.0)
y.append(75.0)
y.append(217.5)
y.append(0.0)
y.append(100.0)
y.append(150.0)
y.append(200.0)
y.append(50.0)
y.append(125.0)
y.append(175.0)
y.append(25.0)
y.append(75.0)
y.append(295.0)
y.append(280.0)
y.append(0.0)
y.append(100.0)
y.append(150.0)
y.append(50.0)
y.append(262.5)
y.append(125.0)
y.append(25.0)
y.append(75.0)
y.append(237.5)
y.append(357.5)
y.append(0.0)
y.append(100.0)
y.append(50.0)
y.append(342.5)
y.append(212.5)
y.append(25.0)
y.append(325.0)
y.append(75.0)
y.append(187.5)
y.append(0.0)
y.append(300.0)
y.append(50.0)
y.append(162.5)
y.append(25.0)
y.append(275.0)
y.append(137.5)
y.append(0.0)
y.append(250.0)

z.append(0.0143170002848)
z.append(0.138046994805)
z.append(0.872992992401)
z.append(0.00615999987349)
z.append(0.154147997499)
z.append(0.828154027462)
z.append(0.0343139991164)
z.append(0.307231009007)
z.append(0.89982098341)
z.append(0.0175800006837)
z.append(0.0719190016389)
z.append(0.478397995234)
z.append(0.93592697382)
z.append(5.99999984843e-05)
z.append(0.0520699992776)
z.append(0.724943995476)
z.append(0.965824007988)
z.append(1.49999996211e-05)
z.append(0.233434006572)
z.append(0.00965300016105)
z.append(0.751917004585)
z.append(0.988830983639)
z.append(0.511997997761)
z.append(0.00044999999227)
z.append(0.00711099989712)
z.append(0.938435018063)
z.append(0.983743011951)
z.append(0.000253000005614)
z.append(0.0598109997809)
z.append(0.50435000658)
z.append(0.00129499996547)
z.append(0.964758992195)
z.append(0.997349977493)
z.append(0.224617004395)
z.append(0.736204981804)
z.append(0.000510999991093)
z.append(0.00981699954718)
z.append(0.937159001827)
z.append(0.00227000005543)
z.append(0.0664549991488)
z.append(0.288078010082)
z.append(0.771333992481)
z.append(0.00428200000897)
z.append(0.11282800138)
z.append(0.319292992353)
z.append(0.00520000001416)
z.append(0.0331760011613)
z.append(0.998526990414)
z.append(0.968023002148)
z.append(0.0103479996324)
z.append(0.0611589998007)
z.append(0.25580099225)
z.append(0.0310849994421)
z.append(0.877821028233)
z.append(0.115939997137)
z.append(0.0475809983909)
z.append(0.0492069981992)
z.append(0.716596007347)
z.append(0.998637974262)
z.append(0.0647180005908)
z.append(0.119029000401)
z.append(0.0855640023947)
z.append(0.98787599802)
z.append(0.457558006048)
z.append(0.11486300081)
z.append(0.925565004349)
z.append(0.14404900372)
z.append(0.339305013418)
z.append(0.2295909971)
z.append(0.844295024872)
z.append(0.246871992946)
z.append(0.346192002296)
z.append(0.294302999973)
z.append(0.655534029007)
z.append(0.360406011343)
z.append(0.347315013409)
z.append(0.621686995029)

h_parent = TH2F(’’,’’,100,min(x)-1,max(x)+1,100,min(y)-1,max(y)+1)
#h_parent = TH2F(’’,’’,100,0.,800.,100,0.,500.)
h_numbers = h_parent.Clone()
for i,element in enumerate(x):
print i,x[i],y[i],z[i]
h_numbers.Fill(x[i],y[i],z[i])
g = TGraph2D(len(x),x,y,z)
g.SetHistogram(h_parent)
h_colors = g.GetHistogram()
[/code]
Above code leads to a segfault in my case (I also tried on another machine, same result).

Some empiric fun facts:

If I remove any (?) element from the arrays, e.g.

x.remove(x[22]) y.remove(y[22]) z.remove(z[22]) it works.

If I add any (?) element to the array, e.g.

x.append(500.) y.append(150.) z.append(.5) it works.

If I do both of the above, i.e.

x.remove(x[22]) y.remove(y[22]) z.remove(z[22]) x.append(500.) y.append(150.) z.append(.5) it does not work.

If I change the borders of the TH2F, i.e.

h_parent = TH2F('','',100,-100.,800.,100,0.,500.) it may work, depending on the values of xlow, xup, ylow, yup.

I fail to make any sense of this. Help is much appreciated.

Thanks,
Basil

Basil,

sorry for the late reply, but I’m traveling (and will be for quite a bit).

Looks to me like a bug in TGraphDelaunay.cxx. I’d argue for:[code]diff --git a/hist/hist/src/TGraphDelaunay.cxx b/hist/hist/src/TGraphDelaunay.cxx
index 8c867c9…3cc2a5c 100644
— a/hist/hist/src/TGraphDelaunay.cxx
+++ b/hist/hist/src/TGraphDelaunay.cxx
@@ -252,7 +252,7 @@ L1:
if (swap) goto L1;

// expand the triangles storage if needed
  • if (fNdt> fTriedSize) {
  • if (fNdt>=fTriedSize) {
    Int_t newN = 2*fTriedSize;
    Int_t *savep = new Int_t [newN];
    Int_t *saven = new Int_t [newN];[/code]
    If you concur, please file a report under http://root.cern.ch/bugs.

Cheers,
Wim

Thanks for the input. I will try to reproduce this using CINT and see if doing >= instead of > fixes the problem.

with CINT I get the Seg Fault when I quit ROOT.
The fix suggested does not help.

The following macro gives the Seg Fault on quit. The faulty line is SetHistogram(). I am investigating.

Hi,

maybe it does not help, but if I read the code, I see allocation: fPTried = new Int_t[fTriedSize]; fNTried = new Int_t[fTriedSize]; fMTried = new Int_t[fTriedSize];
and the use of fNdt:[code] if (fNdt> fTriedSize) {

}

// store a new Delaunay triangle
fNdt++;
fPTried[fNdt-1] = ps;
fNTried[fNdt-1] = ns;
fMTried[fNdt-1] = ms;[/code]
then I still argue for ‘>=’. (In my case, gdb complained about segfault in _int_free in the block, which was gone after the code change). I.e. if fNdt == fTriedSize the block will not be reallocated and the entries will be accessed at fTriedSize ( (fNdt +1 -1) == fTriedSize ), which is one too far.

Note that that same allocation pattern is once more in the block itself.

Thanks for the CINT code. I’ll give it a spin as well (and yes, also in python, the GetHistogram is the problem as it calls TGraphDelaunay::FileIt() eventually).

Cheers,
Wim

I will apply your patch. But that does not fix this issue.

Hi,

well FWIW, in my case I see w/ the CINT code:#8 <signal handler called> #9 0x00007f0fb3d4dff4 in _int_free () from /lib64/libc.so.6 #10 0x00007f0fb1cf4476 in TGraphDelaunay::FileIt(int, int, int) () from /home/wlav/rootdev/v5-34-00-patches/lib/libHist.so #11 0x00007f0fb1cf660e in TGraphDelaunay::Interpolate(double, double) () from /home/wlav/rootdev/v5-34-00-patches/lib/libHist.so #12 0x00007f0fb1cf6d78 in TGraphDelaunay::ComputeZ(double, double) () from /home/wlav/rootdev/v5-34-00-patches/lib/libHist.so #13 0x00007f0fb1cfaaab in TGraph2D::GetHistogram(char const*) () from /home/wlav/rootdev/v5-34-00-patches/lib/libHist.sowhen using ‘>’, whereas it runs fine when using ‘>=’.

Best I can think of …

Cheers,
Wim

In my case, with CINT, running on Mac the macro I sent you, I get no error at run time. I get only an seg fault when I do “.q” and the messages are:

root [0] .x graph2d14.C
Info in <TCanvas::MakeDefCanvas>:  created default TCanvas with name c1
root [1] .q

 *** Break *** segmentation violation



===========================================================
There was a crash.
This is the entire stack trace of all threads:
===========================================================

Thread 1 (process 41067):
#0  0x00007fff871756ac in wait4 ()
#1  0x00007fff8937903a in system ()
#2  0x0000000107c4c08f in TUnixSystem::StackTrace ()
#3  0x0000000107c49ed7 in TUnixSystem::DispatchSignals ()
#4  <signal handler called>
#5  0x000000010962d90d in TGraph2D::Clear ()
#6  0x000000010962d721 in TGraph2D::~TGraph2D ()
#7  0x000000010962d5df in TGraph2D::~TGraph2D ()
#8  0x0000000107bfdcc4 in TCollection::GarbageCollect ()
#9  0x0000000107c01794 in TList::Delete ()
#10 0x0000000108f356c4 in TDirectoryFile::Close ()
#11 0x0000000108f40a02 in TFile::Close ()
#12 0x0000000107bb5c2d in (anonymous namespace)::R__ListSlowClose ()
#13 0x0000000107bb575c in TROOT::CloseFiles ()
#14 0x0000000107c4b8e0 in TUnixSystem::Exit ()
#15 0x0000000107b7796c in TApplication::ProcessLine ()
#16 0x0000000108b70cca in TRint::HandleTermInput ()
#17 0x0000000108b6f367 in TTermInputHandler::Notify ()
#18 0x0000000108b7153d in TTermInputHandler::ReadNotify ()
#19 0x0000000107c493a6 in TUnixSystem::CheckDescriptors ()
#20 0x0000000107c48bf1 in TUnixSystem::DispatchOneEvent ()
#21 0x0000000107bd06ea in TSystem::InnerLoop ()
#22 0x0000000107bd051d in TSystem::Run ()
#23 0x0000000107b789c4 in TApplication::Run ()
#24 0x0000000108b70660 in TRint::Run ()
#25 0x0000000107b6fa5f in main ()
===========================================================


The lines below might hint at the cause of the crash.
If they do not help you then please submit a bug report at
http://root.cern.ch/bugs. Please post the ENTIRE stack trace
from above as an attachment in addition to anything else
that might help us fixing this issue.
===========================================================
#5  0x000000010962d90d in TGraph2D::Clear ()
#6  0x000000010962d721 in TGraph2D::~TGraph2D ()
#7  0x000000010962d5df in TGraph2D::~TGraph2D ()
#8  0x0000000107bfdcc4 in TCollection::GarbageCollect ()
#9  0x0000000107c01794 in TList::Delete ()
#10 0x0000000108f356c4 in TDirectoryFile::Close ()
#11 0x0000000108f40a02 in TFile::Close ()
#12 0x0000000107bb5c2d in (anonymous namespace)::R__ListSlowClose ()
#13 0x0000000107bb575c in TROOT::CloseFiles ()
#14 0x0000000107c4b8e0 in TUnixSystem::Exit ()
===========================================================


Root > 

Using “>” or “>=” in TGraphDelaunay::FileIt makes no difference.

Hi,

there are two differences then: I was running on Linux, and in batch (remote machine). But neither should matter (or at least, I would not know why they’d matter).

I just had to compile ROOT on my Mac laptop for another topic on this forum, so I gave it a try: it works perfectly fine either way! :slight_smile:

However, I don’t see a default TCanvas being created? (No message nor TCanvas popping up.)

Cheers,
Wim

Oops, yes sorry, in my local version I added a g->Draw() to the macro I posted here (i like to see drawing :slight_smile: )
But that doesn’t make any difference, I get the Seg Fault on .q with or without this Draw().

I added a protection in TGraph2D::Clear() to avoid this Seg Fault on .q. Now in 5.34 and trunk.

Hi

Thanks for investigating. The above macro works for me, but if I add a g->Draw();it crashes. The crash also persists on the trunk.

Basil

I do not get any crash with the trunk and with 5.34 using a C macro. I am running on Mac. I added the Draw() as you suggested.
graph2d14.C (2.74 KB)

Same macro but with crash on

However, if I apply the fix suggested by Wim, it runs as it should.

[code]diff --git a/hist/hist/src/TGraphDelaunay.cxx b/hist/hist/src/TGraphDelaunay.cxx
index 8c867c9…3cc2a5c 100644
— a/hist/hist/src/TGraphDelaunay.cxx
+++ b/hist/hist/src/TGraphDelaunay.cxx
@@ -252,7 +252,7 @@ L1:
if (swap) goto L1;

// expand the triangles storage if needed
  • if (fNdt> fTriedSize) {
  • if (fNdt>=fTriedSize) {
    Int_t newN = 2*fTriedSize;
    Int_t *savep = new Int_t [newN];
    Int_t *saven = new Int_t [newN];[/code]

Cheers,
Basil

Wim’s fix is now in 5.34 and trunk.