From: Alexandr Malusek <Alexandr.Malusek_at_imv.liu.se>

Date: Tue, 26 Apr 2005 15:58:33 +0200

Date: Tue, 26 Apr 2005 15:58:33 +0200

Hi,

There seems to be a bug in the TGraph::Eval() routine:

Test code:

{

const Int_t dim = 5; Double_t x[dim] = {1.0, 2.0, 3.0, 4.0, 5.0}; Double_t y[dim] = {1.0, 2.0, 3.0, 4.0, 5.0};TGraph *g = new TGraph(dim, x, y);

cout << "y(4.99) = " << g->Eval(4.99) << endl; cout << "y(5.00) = " << g->Eval(5.00) << endl; cout << "y(5.01) = " << g->Eval(5.01) << endl;}

Output:

y(4.99) = 4.99 # correct y(5.00) = 1.01246e+07 # wrong y(5.01) = 1.01449e+07 # wrong

I think the code in root/graf/src/TGraph.cxx:

Double_t TGraph::Eval(Double_t x, TSpline *spline, Option_t *option) const

...

* //linear interpolation*

* //find point in graph immediatly below x*

* //In case x is < fX[0] or > fX[fNpoints-1] return the extrapolated point*

Int_t up, low = TMath::BinarySearch(fNpoints,fX,x); if (low == fNpoints-1) {up=low; low = up+1;} if (low == -1) {low=0; up=1;} if (fX[low] == fX[low+1]) return fY[low]; Double_t yn = x*(fY[low]-fY[low+1]) +fX[low]*fY[low+1] - fX[low+1]*fY[low]; return yn/(fX[low]-fX[low+1]);

should be

* //linear interpolation*

* //find point in graph immediately below x*

* //In case x is < fX[0] or > fX[fNpoints-1] return the extrapolated point*

Int_t low = TMath::BinarySearch(fNpoints,fX,x); Int_t up = low+1; if (low == fNpoints-1) {up=low; low = up-1;} if (low == -1) {low=0; up=1;} if (fX[low] == fX[up]) return fY[low]; Double_t yn = x*(fY[low]-fY[up]) +fX[low]*fY[up] - fX[up]*fY[low]; return yn/(fX[low]-fX[up]);

The problem is with

if (low == fNpoints-1) {up=low; low = up+1;}

the rest are cosmetic changes to utilize the variable "up". The modified code gives correct results:

y(4.99) = 4.990000e+00 y(5.00) = 5.000000e+00 y(5.01) = 5.010000e+00

The problem was observed with root-4.03/02 and root-4.03/04 on Linux, other versions may also be affected.

Regards,

Alexandr
Received on Tue Apr 26 2005 - 15:58:41 MEST

*
This archive was generated by hypermail 2.2.0
: Tue Jan 02 2007 - 14:45:07 MET
*