How do I handle <Double_t buffer> objects?

Hello PyRoot folks,

I have nearly completed the conversion of one of my original ROOT scripts into python, except for one impasse. Within my fit routine, I have having difficulty translating the following code into python.

Double_t fitFunction(Double_t *x, Double_t *par){   
    Double_t value = quadBackground(x, &par[3]) + gaussPeak(x, par);
    return value;
}

I learned that par is a buffer object instead of a list like I originally thought, so slices won’t work. I’ve also messed around with the built-in buffer module, but without success.

I’m sure many of you have had this problem and have found an easy solution. If so, would you mind sharing it?

Thank you,

Sean

1 Like

Sean,

assuming that I understand your situation correctly, I can’t think of anything better than:a = array.array( 'd', [ par[3] ] ) value = quadBackground(x, a) + gaussPeak(x, par); par[3] = a[0] return value
Which version of ROOT or you using?

Cheers,
Wim

I’m using the latest version of ROOT and python, 5.16 and 2.5.

I hashed together a quick sample that gives a similar error.

from ROOT import TMath, TCanvas, TFile, TF1
import array

def quadratic_background_gaussian_peak(x, par):
   #I wished this worked!
   par2 = array.array('d', [par[3]] ) # Workaround starts 
   value =  single_gauss(x, par) + polynomial_order_2(x, par2)
   par[3] = par2[0]
   return value
   # I don't want to do this.
   #return par[0]/(TMath.Sqrt(2*TMath.Pi())*par[2]) * TMath.Exp( -0.5*((x[0]-par[1])/par[2])**2 ) + \
   #     par[3] + par[4]*x[0] + par[5]*x[0]*x[0]
   
def single_gauss(x, par):
   return par[0]/(TMath.Sqrt(2*TMath.Pi())*par[2]) * TMath.Exp( -0.5*((x[0]-par[1])/par[2])**2 )
   
def polynomial_order_2(x, par):
   return par[0] + par[1]*x[0] + par[2]*x[0]*x[0]

mycanvas = TCanvas("c1","mycanvas",800,600)
myfile = TFile("spec_dop.root")
myhisto = myfile.Get("h1")
myhisto.Rebin(10)
myfit = TF1("myfit",quadratic_background_gaussian_peak,995.,1220.,6)
myfit.SetParameters(40, 1112.0, 10.0, 4.0, 0.0, 0.0, 0.0)
myhisto.Fit(myfit,"R+")

I’ve sent the very small data file as an attachment.

The simple solution would be to write out everything, but in a few of my classes I have written functions to dynamically generate TF1 objects, and it would be a pain to write every possible function combination I might use!

Sean
spec_dop.root (6.39 KB)

Hi,

just a quick reply, as I appear to have misunderstood (I’ll try your code tomorrow and look in more detail): if you only want to slice the array, and don’t care about modifying it or passing it to a C++ function (which is what I originally thought), just copying it into a list and slicing that:slice = list(par)[2:]will do the trick.

The slicing of the buffer directly fails b/c it is sliced on 1-byte boundaries instead of sizeof(double) boundaries and the returned as a string. I’ve put fixing that on the TODO list (ATLAS software workshop, next week, so busy :confused: ).

Cheers,
Wim

Hi Wil,

You suggestion worked perfectly.

Thanks!

Sean

Hi,

a related question I believe, maybe simply the same?

I am working in pyroot, and have a TGraph. From the TGraph I would like to get hold of the arrays that make the x and y-axis, and feed these into a new TGraphAsymmError object. I tried the following:

x_arr = graph.GetX()
y_arr = graph.GetY()

This seems to be an object of type buffer.

I want to use these arrays again to pass into a TGraphAsymmError object. How to convert the arrays properly?

I tried

x_list = list(x_arr)
x_new_array = asarray(x_list)

But I get a segmentation violation at list(x_arr), and I also find this “solution” a bit clumsy.

There must be a straightforward way?

Thanks!
Maiken

Maiken,

the returns from GetX() and GetY() are buffers b/c the C++ header only declares a Double_t*, from which no further information (such as the size of the array) can be had. If you want to pass that buffer to a python list, you’ll first need to fix up its size (I think that the call GetN() on the TGraph hands the right size; use SetSize(newsize) on the buffer to set its size). The reason being that the list iterates over the buffer to copy the values, and so w/o its size being set, it will iterate beyond the end of the underlying array.

However, I’m puzzled as to why you want to construct a new array to pass into TGraphAsymmError? AFAICS, it takes a Double_t* in its ctor, so the buffer can be readily passed in, no?

Cheers,
Wim

Thank you, excellent. I did not understand that it could just be used like that.

Works very well, thank you.

I wonder about a second thing: Really what I would like is to port the values from a histogram into the TGraphAsymmErrors, but I did not find a way to elegantly get both the x- and y-values. So I first made the histogram into a TGraph and then used the GetX() and GetY().

Am I right that one cannot grab the array of x and y-values directly from the histogram itself?

Thanks a lot,
Maiken

Maiken,

the y-values can be hadd from the GetArray() function that histograms inherited from TArrayF/D. The x-values I think would have to be build up from a new array.

(I see that TGraphAsymErrors has a contructor that takes a TH1*, but it loops over the errors of the histogram, not the values, so I presume that’s not what you want.)

Cheers,
Wim