Hello, I measured some physics quantities values using 3 different methods then I have the three values, for example one is the ratio of two parameter and I got :
a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036
is there an automatic way to calculate the according between the three result (For example how many sigma)?
Thank you
Please read tips for efficient and successful posting and posting code
ROOT Version: Not Provided
Platform: Not Provided
Compiler: Not Provided
Try searching the forum for that. Here is one example result
Thank you @dastudillo I read @moneta 's code, but I didn’t understand so much…
std::vector<double> x(100);
gRandom->RndmArray(x.size(),x.data());
TMath::StdDev(x.begin(), x.end() );
for example, to calcucalte how many sigma is the according between the three results
a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036
where do I have to write these number in @moneta 's code?
I tried to write
gRandom->RndmArray(x.size(3),x.data(0.692, 0.6833, 0.6712));
but I got this error
I also tried @eguiraud’s code in this way
using namespace ROOT::VecOps;
RVec<float> rv = {0.11, 0.0053, 0.0036};
auto srddev = StdDev(rv);
but
a. Does this methos use only the errors of the results and not the values (0.692, 0.6833, 0.6712)?
b. You see…ROOT closes automatically when I run the
auto srddev = StdDev(rv);
thank you
Yes, RVec::StdDev does a “dumb” calculation of the square root of the variance of the population, nothing fancy
That really should not happen. In fact I typed those exact same commands in a ROOT prompt on my laptop and things are fine:
~ root -l
root [0] using namespace ROOT::VecOps;
root [1] RVec<float> rv = {0.11, 0.0053, 0.0036};
root [2] auto srddev = StdDev(rv);
root [3] srddev
(double) 0.060945248
root [4]
Not sure what’s going on there 
Thank you @eguiraud I tried the code using gentoo and it worked
I use Windows… @bellenot could it be a bug on Root for windows?
@eguiraud you said that this way is a dumb way… then, do you think that it’s better @moneta’s code
std::vector<double> x(100);
gRandom->RndmArray(x.size(),x.data());
TMath::StdDev(x.begin(), x.end() );
to get the according between the 3 measurements?
@moneta please, can you check the problem? I resume… I measured physics quantities using 3 instruments…for example for a ratio I got these values (one for each instrument)
a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036
I need a test to estimate the according between the 3 results.
I followed your code
std::vector<double> x(100);
gRandom->RndmArray(x.size(),x.data());
TMath::StdDev(x.begin(), x.end() );
then I wrote
std::vector<double> x(100);
gRandom->RndmArray(x.size(3),x.data(0.692, 0.6833, 0.6712));
but I got this error
Can you check please what I did wrong?
and, is there a better way for my goal how many sigma is the according?
Thank you
No no I meant that RVec::StdDev is dumb, i.e. it does not do anything smart like taking into account error values attached to the inputs 
As far as I understand TMath::StdDev does the same by default but it also has an overload that can take into account input weights, you might be able to do something with that.
This is not a valid syntax (just try on Linux)
And this code works for me on Windows:
C:\Users\bellenot\DAQ>root -l
root [0] using namespace ROOT::VecOps;
root [1] RVec<float> rv = {0.11, 0.0053, 0.0036};
root [2] auto srddev = StdDev(rv);
root [3] srddev
(double) 0.060945248
root [4]
@eguiraud yes, I understood that you meant the RVec::StdDev is dumb, that’s because I asked to @moneta if there is a better way for my goal!
@moneta: important… I said that I made 3 measurements using three instruments…but I got results using the measurements and making some calucations and fit …
then, for example
a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036
are NOT the values read on the instruments, but are results of calculations using the values read by the instruments
I don’t know, try to follow what @eguiraud or @dastudillo proposed…
Thank you @bellenot I tried to follow the example showed by @dastudillo… but @eguiraud’s code doesn’t look like what I need… I was trying @moneta’s code but you see, I dont’ understand how to apply it…so I hope @moneta can help me, I know he is the most suitable to help regarding on the analysis way using ROOT…
or maybe, also @jblomer could know a way to calculcate the according between the results (because of he works in particle physics)
Hi,
I am not sure I have understood you completly.
You have 3 measurements, and what you would like to know, if they agree to each other ?
If this is the case, it is very simple, you can compute the mean, and compute their chisquared value, from which you can obtained a p-value. To do this, is easier by performing a simple fit with an horizontal line:
{
std::vector<double> x = {1,2,3};
std::vector<double> y = { 0.692, 0.6833, 0.6712 };
std::vector<double> ey = { 0.011, 0.0053, 0.0036};
auto g = new TGraphErrors(3,x.data(),y.data(), nullptr, ey.data());
auto result = g->Fit("pol0","S");
double chi2 = result->Chi2();
double ndf = result->Ndf();
double pvalue = ROOT::Math::chisquared_cdf_c(chi2, ndf);
std::cout << "p value of the fit : " << pvalue << std::endl;
// number of standard sigma. Divide p-value by 2 since deviations can be in both sides
std::cout << "number of sigma is " << ROOT::Math::normal_quantile_c(0.5*pvalue,1) << std::endl;
}
Lorenzo
Thank you @moneta it was what I needed!
@bellenot I tried the code on my root (windows machine) but when I lunch the command
double ndf = result->Ndf();
Root crashes
I also tried the macro
sigmavalue.cpp (2.6 KB)
but I get this error
Is it a problem on root for windows?
Thank you
It works for me with ROOT master (6.25/01):
C:\Users\bellenot\rootdev>root
------------------------------------------------------------------
| Welcome to ROOT 6.25/01 https://root.cern |
| (c) 1995-2021, The ROOT Team; conception: R. Brun, F. Rademakers |
| Built for win32 on Mar 09 2021, 17:28:03 |
| From heads/master@v6-25-01-74-g202a84bb2e |
| With MSVC 19.24.28319.0 |
| Try '.help', '.demo', '.license', '.credits', '.quit'/'.q' |
------------------------------------------------------------------
root [0] std::vector<double> x = {1,2,3};
root [1] std::vector<double> y = { 0.692, 0.6833, 0.6712 };
root [2] std::vector<double> ey = { 0.011, 0.0053, 0.0036};
root [3] auto g = new TGraphErrors(3,x.data(),y.data(), nullptr, ey.data());
root [4] auto result = g->Fit("pol0","S");
****************************************
Minimizer is Linear / Migrad
Chi2 = 5.78669
NDf = 2
p0 = 0.67618 +/- 0.0028745
root [5] double chi2 = result->Chi2();
root [6] double ndf = result->Ndf();
root [7] double pvalue = ROOT::Math::chisquared_cdf_c(chi2, ndf);
root [8] std::cout << "p value of the fit : " << pvalue << std::endl;
p value of the fit : 0.0553906
root [9] std::cout << "number of sigma is " << ROOT::Math::normal_quantile_c(0.5*pvalue,1) << std::endl;
number of sigma is 1.9158
root [10]
And if you look at the (many) compilation errors with your code, you should be able to fix them. If you want me to try, please provide also the input file.
You are at least missing several headers:
#include "TF1.h"
#include "Math/ProbFuncMathCore.h"
#include "Math/QuantFuncMathCore.h"
#include <cstdlib>
#include <iostream>
#include <fstream>
Then it compiles fine:
C:\Users\bellenot\rootdev>root -l
root [0] .L sigmavalue.cpp+
Info in <TWinNTSystem::ACLiC>: creating shared library C:/Users/bellenot/rootdev/sigmavalue_cpp.dll
sigmavalue_cpp_ACLiC_dict.cxx
Creating library C:/Users/bellenot/rootdev\sigmavalue_cpp.lib and object C:/Users/bellenot/rootdev\sigmavalue_cpp.exp
sigmavalue_cpp_ACLiC_dict.cxx
Creating library C:/Users/bellenot/rootdev\sigmavalue_cpp.lib and object C:/Users/bellenot/rootdev\sigmavalue_cpp.exp
root [1] .q
Thank you @bellenot adding the headers the macro worked on windows too!
@moneta sorry…one more question please…
to quantify the agree between an experimental measurement and the theoretical value, is there some function on ROOT?
For example, the Planck’s constant value is h=6.62607004 *10^-34Js…
I got h=(6.424+/- 0.052)*10^-34 Js
I want to quantify the goodness of the experimental result… is there a function in ROOT?
Thank you
Hi,
In this case with a single measurement, you just compute the difference divide by the error.
This will give you number of sigma’s.
Lorenzo
Thank you @moneta
This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.